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Spectral  Analysis  via  Quadratic  Frequency-Smoothing 
of  Fourier-Transformed,  Overlapped,  Weighted  Data  Segments 

Introduction 

Spectral  analysis  techniques  have  received  a  great  deal  of  attention  in  the  past 
(references  1-12),  ranging  from  the  original  autocorrelation  approach  of  Blackman- 
Tukey  (reference  2)  to  the  more  recent  weighted,  overlapped,  segment-averaging 
FFT  approach  (references  7-12).  These  two  apparently  disparate  approaches  are 
shown  here  to  be  limiting  special  cases  of  a  generalized  framework  for  spectral 
analysis;  thus  consideration  of  this  general  technique  elucidates  the  fundamental 
behavior  and  performance  of  a  rather  wide  variety  of  spectral  approaches  and  their 
tradeoffs.  This  generalized  framework  has  already  been  presented  in  references  13- 
15,  where  a  brief  summary  of  some  of  the  main  features  has  been  mentioned. 
Additionally,  some  of  the  analytical  results  to  be  presented  here  were  alluded  to 
there;  however,  none  of  the  detailed  derivations  or  quantitative  results  in  this  report 
were  given  at  that  time. 

There  are  two  fundamental  parameters  that  critically  affect  the  performance  of 
any  spectral  estimation  technique.  They  are  the  available  record  length,  T,  of  the 
stationary  random  process  under  investigation,  and  the  effective  frequency 
resolution,  B^,  of  the  technique  under  consideration.  We  would  like  to  be  able  to 
attain  fine  resolution  (small  Bc)  with  short  data  lengths  and  storage  (small  T); 
however,  stable  results  (small  fluctuations)  are  achievable  only  if  the  product  TBC  is 
much  larger  than  unity.  The  problems  we  address  are  how  to  make  optimum  use  of 
a  given  limited  amount  of  data  in  order  to  realize  a  specified  desired  resolution  w  ith 
maximum  stability,  and  to  determine  what  tradeoffs  are  available  regarding  win¬ 
dowing  and  weighting  at  different  stages  of  the  spectral  analysis  procedure.  It  is 
assumed  that  the  reader  is  familiar  with  the  tradeoffs  presented  in  reference  9  for 
the  weighted,  overlapped,  segment-averaging  FFT  procedure. 

The  generalized  framework  for  spectral  analysis  that  is  presented  here  is  capable 
of  a  wide  variety  of  forms  in  addition  to  the  Blackman-Tukey  and  FFT  approaches 
mentioned  above.  In  order  to  compare  these  various  forms  with  each  other  on  a 
reasonable  basis,  it  is  required  that  each  analysis  technique  realize  the  same  effective 
resolution  bandwidth,  Be,  and  that  they  all  utilize  the  same  data  record  length  T. 
Without  these  reasonable  constraints,  valid  conclusions  about  relative  per¬ 
formances  of  different  techniques  are  tentative  at  best.  This  insistence  upon  equal 
effective  frequency  resolution  necessitates  a  rather  detailed  investigation  of  the 
effects  of  the  weightings  and  windows  employed  in  the  generalized  framework  and 
their  allowed  durations.  The  desirability  of  an  overall  effective  window  for  spectral 
analysis  with  low  side  lobes  and  good  decay  is  achievable  only  through  careful 
choice  of  the  combined  weightings.  The  constraint  upon  the  effective  frequency 
resolution  naturally  also  shows  up  in  the  analysis  of  the  variance  of  the  spectral 
estimation  technique,  as  well  as  in  its  average  value,  leading  to  some  numerical 
analysis  complications;  nevertheless,  it  is  believed  to  be  the  proper  basis  of  com¬ 
parison  and  is  maintained  throughout. 
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The  two  major  statistical  parts  of  this  report  deal  with  the  mean  and  the  variance 
of  the  spectral  estimate.  The  result  for  the  average  value  leads  to  the  definition  of 
the  effective  window  of  the  generalized  spectral  analysis  technique,  in  terms  of  the 
temporal  and  lag  windows.  The  variance  result  incorporates,  additionally,  the 
amount  of  overlap,  the  number  of  data  pieces,  and  the  ambiguity  functions  of  the 
temporal  and  lag  windows;  the  complexity  of  the  latter  results  debilitates  easy  in¬ 
terpretation  and  it  has  been  found  necessary  to  resort  to  numerical  evaluation  of  the 
variance,  for  practical  cases  of  interest. 


2 


Ultimate  Stability  Attainable  From  a  Given  Record  Length 

Suppose  a  stationary  (complex)  data  record  x(t)  of  length  T  seconds  is  available, 
and  that  we  wish  to  estimate  its  power  density  spectrum*  G(f)  with  an  effective 
frequency  resolution  of  Bc  Hz,  where  WQ(f)  is  the  narrowband  window  through 
which  the  power  density  spectrum  is  to  be  observed.  These  two  frequency-domain 
quantities  are  related  according  tot 

B  .  IidfA(fji 

'  /«  "o<«  (1) 

This  bandwidth  measure,  Be,  is  called  the  statistical  bandwidth  of  W0(f)  in  reference 
5,  page  265.  The  relation  of  effective  bandwidth  Be  to  half-power  bandwidth  Bh  is 
considered  in  appendix  A;  it  is  shown  that  for  good  windows,  the  ratio  of  the  two 
bandwidths  is  relatively  independent  of  the  exact  window  shape.  Thus  it  is  possible 
to  translate  results  to  other  bandwidth  measures  without  significantly  affecting  the 
essential  quantitative  aspects. 

If  we  take  the  original  data  record  and  pass  it  through  a  narrowband  linear 
(complex)  filter  with  power  transfer  function  equal  to  the  window,  |H(f)|2  =  Wo(f), 
and  which  is  centered  at  a  frequency,  fQ,  of  interest,  we  will  have  lost  no  relevant 
information  about  the  process  in  the  frequency  band  of  interest,  because  we  have 
filtered  out  information  of  no  use.  We  can  now  estimate  the  power  in  the 
narrowband  filter  output  process  and  use  it  as  a  measure  of  the  spectrum  of  the 
input  process  in  the  neighborhood  of  frequency  f0.  See  figure  1 . 


Let  z(t)  be  the  complex  output  process  from  the  narrowband  filter  when  excited 
by  the  available  T  seconds  of  data  x(t).  If  we  ignore  a  starting  transient  (i.e. , 
assuming  T  »  1/Be),  the  filter  output  power  estimate  in  the  band  of  w  idth  Be  is 

P  =  Y  fdt  J  z(t) | 2  =  /  dt  g(t)  |z(t)|2  ,  (2) 


•For  brevity,  we  use  the  term  spectrum  rather  than  autospectrum  in  this  report, 
flntegrals  without  limits  are  over  the  range  of  the  nonzero  integrands. 
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where  gate  function 


g(t) 


1/T  for  t  e  T 
0  otherwise 


(3) 


The  measure  of  stability  we  adopt  for  this  estimator,  and  for  the  others  to  follow, 
is  the  quality  ratio  defined  as  2 


Q  2 


A 

Var(P) 
Av2 (P) 


£2 


A 

P 


A. 

P 


(4) 


where  Av(£)  and  Var(P)  denote  the  average  value  and  variance  of  P,  respectively, 
and  an  overbar  denotes  an  ensemble  average.  We  have  average  value 

Av(P)  -  /dtg(t)  |z(t)  |2  =  |z(t)|2  =  Rz(0)=  Jdf  G(f)  |H(f )  | 2 


=  G(fo)  / df  |H(f)|2  ,  (5) 

assuming  that  filter-input  spectrum  G(f)  does  not  vary  quickly  with  respect  to  Be,  in 
the  neighborhood  of  f0.  Rz(x)  is  the  correlation*  of  filter  output  process  z(t). 

Also,  we  have  mean  square  value 

P2  =  J  J  dt  du  g(t)  g(u)  |z(t)|2  |  z  (u)  | 2  .  (6) 


Now  in  the  interval  T,  filter  output  z(t)  will  be  approximately  a  stationary  zero- 
mean,  complex,  analytic  Gaussian  process  for  small  Be;  filter  H(f)  has  filtered  out 
zero  and  all  negative  frequencies.  Then  fourth-order  moment 


z(t)z*(t)z(u)z*(u)  ?R2(0)  +  |R  (t  -  u)|2  . 

*  P) 

There  follows  from  (4)  and  (6), 

Var(P)  =  f  J  dt  du  g(t)g(u)  | Rz Ct  -  u)  | 2  =  J dx  <t>g(x)  |Rz(x)  |2  , 

(8) 


where  gate-correlation  t  of  function  g(t)  is 

<J>  (r)  =  f dt  g(t)g(t  -  t) 

s  •* 


(9) 


•For  brevity,  we  use  the  term  correlation  instead  of  autocorrelation  in  this  report. 

tFor  stationary  processes,  we  let  R  denote  the  ensemble-average  correlation,  whereas  for  aperiodic 

nonrandom  functions,  we  let  $  denote  the  integral  correlation;  see  (5)  and  (9). 
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Since  the  gate-correlation  +g(T)  extends  over  ±T,  while  process  correlation  Rz(t) 
extends  only  over  approximately  ±3/Be,  we  have,  via  Parseval’s  theorem  and  for 
TBe»  1, 

Var(P)  s  4>g(0)  /dx|R2(T]|2  =  I  f  df  G^ff) 

=  I  f  df  |H(f)|4  G2(£)  *  I  G2(f  )  fdf  |H(f)  | 4  . 

1  J  i  o  j  (JO) 

The  quality-ratio  measure  of  stability  is  then,  from  (4),  (5),  ( 1 0),  and  ( 1 ), 
/df|H(f)|4  J df  W2(f) 

Q  =  ir7 - :tt  s  t  --  -  li  -T r  • 

|/<Jf|H(f)  |  J  [/«  "  (||) 

This  is  the  limiting  (smallest)  value  of  Q  for  specified  frequency  resolution  Be  and 
available  record  length  T  when  TBe »  1 .  No  other  spectral  procedures  can  improve 
on  it;  they  can  merely  approximate  it.  As  such,  (1 1)  is  the  benchmark  against  which 
other  procedures  can  be  compared,  under  the  condition  that  T  and  Be  are  equal  to 
those  values  for  the  various  procedures  under  consideration. 

The  normalized  quality  ratio  is  defined  as  Q  •  TBC.  Thus  the  normalized  quality 
ratio  can  never  be  smaller  than  unity,  which  value  can  only  be  approached  for  large 
TBe  through  proper  processing  techniques. 
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Description  of  Spectral  Analysis  Technique 


We  begin  by  defining  a  temporal  weighting  function  w,(t)  of  finite  duration  L,; 
that  is. 


w^Ct)  0  only  for  |t|  <  L^/2 


(12) 


As  shown  in  figure  2,  temporal  weighting  w,(t)  is  real,  even,  and  peaked  at  the 
origin.  Although  this  presentation  is  couched  in  terms  of  continuous  functions,  we 
shall  show  shortly  that  it  includes  discrete  digital  processing  as  a  special  case. 


Figure  2.  Temporal  Weighting  w,(t) 


The  available  data  record  is  x(t)  for  0  <  t  <  T;  this  (complex)  random  process  is 
presumed  second-order  stationary  in  that  observation  interval.  We  shift  the  tem¬ 
poral  weighting  by  L,/2  +  pS  and  multiply  it  by  x(t)  to  generate  the  p-th  piece  of 
weighted  data: 

y  (t)  =  x(t)  w  (t  -  -  ps)  for  0  <  p  <  f  -  l  . 

P  1\  2  /  “  ~  (13) 


Here  p  is  an  integer;  if  shift  S  <  L,,  then  yp(t)  and  yp  +  ,(t)  will  overlap  on  the  t-axis. 

The  first-stage  power  density  spectral  estimate  at  frequency  f  is  obtained  by 
averaging  the  magnitude-squared  value  of  the  Fourier  transform  of  data  piece  yp(t), 
over  a  total  of  P  pieces: 

P-1  I  2 

G,  (f)  =  ^  V  f  dt  exp(-i2irft)  y  (t)  for  any  f  . 

1  Pp-0|J  P  04) 

This  procedure  is  the  same  as  that  considered  in  reference  7  and  in  reference  9,  eqs. 
(2)  and  (3).  Since  x(t)  is  available  only  for  0  <  t  <  T,  we  prevent  the  weighting  in  (13) 
from  extending  beyond  that  interval;  mathematically  this  means  that  we  must  have 


L  +  (P  -  1)  S  s.  T 


(15) 
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An  alternative  interpretation  of  (14)  is  very  illuminating.  We  define  the  inverse 
Fourier  transform  of  (14)  as  the  first-stage  (auto)correlation  estimate;  there  follows 
immediately  at  delay  t, 

Rj(t)  =  J  df  exp(i2irf t)  SjCf) 


=  }  £  Jdt  yp(v  - T>  for  a11  T  *  (i6) 

p=0  r 

where  we  have  allowed  random  processes  x(t)  and  yp(t)  to  be  complex.  This  is 
recognized  as  the  average  of  the  sample  correlations  that  can  be  formed  at  delay  t, 
from  each  of  the  P  pieces  of  weighted  data  in  (13).  Since  temporal  weighting  w,(t)  is 
zero  for  |t|  >  L,/2  according  to  (12),  we  see  from  (13)  and  (16)  that 


Rjfx)  =  0  for  |  x  |  >  Lj 


(17) 


The  parameter,  r,  is  called  the  lag  domain  variable,  because  of  the  way  it  appears  as 
a  delayed  time  in  (16).  Equation  (16)  (and  those  to  follow)  is  true  for  all  t.  Both 
sides  of  (16)  are  zero  over  most  of  the  range  of  t;  nevertheless,  it  is  mathematically 
convenient  to  employ  the  equality  of  both  sides  of  (16)  for  all  r  in  various  trans¬ 
formations  below. 

The  second-stage  power  density  spectral  estimate  is  defined  as  a  frequency- 
smoothed  version  of  the  first-stage  result: 

G2(f)  =  GjCf)  ®  W2(f)  =  / du  G1(u)W2(f  -  u)  , 

where  ®  denotes  convolution.  This  is  termed  quadratic  smoothing  since  it  is  done  in 
terms  of  power  quantities  rather  than  voltages.  Equation  (18)  is  the  desired  output 
from  the  generalized  spectral  analysis  technique  considered  here.  W,(f)  is  called  the 
lag  window,  for  reasons  to  be  given  below.  The  equivalent  statement  to  (18)  in  the 
lag  domain  is  obtained  by  Fourier  transforming  (18);  the  second-stage  correlation  is 


R2(t)  =  J  df  exp(i2irfr)  G2(f)  =  R^r)  w2(t)  ,  ^ 

where  we  used  (18)  and  (16)  and  defined  the  Fourier  transform  pair 
w2(t)  =  J  df  exp(i2wfr)  W2(f)  , 

W2(f)  =  J  dr  exp(-i2Trfx)  w2(r)  .  (20) 

w2(r)  and  W2(f)  are  both  real,  even,  and  peaked  at  their  origins.  Since  w2(r)  appears 
multiplicatively  in  (19),  it  is  called  the  lag  weighting;  its  transform  W,(f)  is  the  lag 
window.  The  convention  adopted  throughout  this  report  is  that  multiplication  by  a 
function  in  the  t  or  t  domains  is  called  a  weighting;  the  counterpart  to  this  operation 
in  the  Fourier  transform  domain  (frequency  f  domain)  is  convolution  and  is  called 
windowing. 
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We  shall  let  lag  weighting  w2(x)  be  of  duration  2L2;  that  is, 

w2(x)  =0  for  | t |  >  L2  .  (21) 

A  typical  plot  is  shown  in  figure  3;  the  reason  for  the  apparent  notational 
discrepancy  between  the  lengths  in  figures  2  and  3  will  become  clear  when  the  lag- 
domain  counterpart  of  temporal  weighting  wt(t)  is  encountered  later. 


Figure  3.  Lag  Weighting  w2(x) 


1 


We  have  already  observed  that  ft  ,(x)  is  zero  for  |x|>L,  in  (17).  Therefore,  it 
follows  from  (19)  and  (21)  that 

R2(t)  =  0  for  | x |  >  min  (Lj,  L2)  .  ^2) 

However,  although  we  must  have  temporal  length  L,  <  T  (from  (15)  for  P  =  1),  there 
is  no  restriction  on  L2.  We  could  have  L2  larger  than  L,  and  T;  this  would  simply 
mean  that  we  would  be  lag-weighting  some  zero  estimates  of  ft,(x)  in  (19)  for  the 
larger  values  of  |x|.  Also  there  are  no  constraints  such  as  realizeability  on  the  lag 
weighting  or  window. 

For  example,  the  special  case  of  no  quadratic  frequency-smoothing  corresponds 
to 

W2(f)  =  6(f),  w2(x)  =  1,  L2  =  <*>  for  no  smoothing, 

(23) 


for  which  (18)  yields  62(f)  =  6,(0-  Thus  we  have  our  standard  first-stage  spectral 
estimate  (14)  as  a  limiting  case  of  the  generalized  spectral  analysis  technique.  On  the 
other  hand,  if  lag  window  W2(f)  were  broad  (small  L2),  there  would  be  a  significant 
amount  of  smoothing  taking  place  in  the  band  about  u  =  f  in  (18)  where  window  W2 
is  non  zero. 

There  is  no  inherent  limitation  on  the  relative  sizes  of  L,  and  L2  as  yet:  L2  can  be 
chosen  as  large  as  desired,  while  L,  is  subject  to  the  upper  bound  T.  However,  when 
we  specify  the  overall  effective  frequency  resolution  of  the  generalized  technique,  a 
relation  between  L,  and  L2  will  ensue. 
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Another  important  special  case  of  the  generalized  spectral  analysis  technique  is 
afforded  by  P  =  1,  w,(t)=  1  for  |t|  <  L,/2,  and  L,  =  T.  Then  (16)  and  (13)  indicate 
that  ft  ,(t)  is  simply  the  sample  correlation  of  the  available  data  x(t)  of  length  T, 
while  ft2(-r)  >n  (19)  is  a  weighted  version  of  ft,(T)  for  |t|  <  L2.  But  this  is  precisely  the 
Blackman-Tukey  approach  described  in  reference  2;  the  choice  of  lag  weighting 
w2(t)  and  its  length  is  fully  discussed  there.  For  example,  if  w2(t)  =  (T— |t|)_i  for 
|t|  <  L2  <  T,  then  ft2(-r)  is  an  unbiased  estimator  for  |r|  <  L2;  see  reference  2,  page 
11. 

For  P>1  and  general  temporal  weighting  w,(t),  lag  weighting  w2(t),  and  overlap, 
a  wide  variety  of  processors  is  possible  via  the  generalized  framework  set  up  above. 
How  should  the  two  weightings  be  traded  off  against  each  other?  Can  the 
deleterious  effects  of  a  poor  or  preselected  temporal  weighting  be  undone  by  proper 
choice  of  lag  weighting?  Recall  that  none  of  these  techniques  can  hope  to  better  the 
quality-ratio  result  (1 1),  but  hopefully,  some  can  do  as  well,  with  less  computational 
effort  and  storage. 

A  related  procedure  to  the  one  presented  here  has  been  given  in  references  16  and 
17.  However,  neither  incorporate  overlapping,  and  the  fundamental  tradeoffs 
between  the  temporal  and  lag  weightings  were  not  studied.  Furthermore,  the  only 
frequency-smoothing  case  considered  was  a  rectangular  boxcar,  which  severely 
limits  the  potential  of  the  technique;  some  advantages  of  the  generalized  technique 
considered  here  will  become  apparent  at  a  later  stage.  For  the  time  being,  we  ob¬ 
serve  that  side  lobe  control  will  be  realized  by  a  mixture  of  temporal  weighting  and 
lag  weighting  (frequency  smoothing),  while  stability  will  be  achieved  by  a  com¬ 
bination  of  segment  averaging  and  frequency  smoothing  (lag  weighting). 

Discrete-Time  Processing 

AH  the  functions  above  have  been  tacitly  assumed  no  worse  than  discontinuous; 
see  figures  2  and  3  for  example.  However,  there  is  nothing  in  the  above  mathematics 
which  precludes  impulsive  behavior.  For  example,  suppose  the  temporal  weighting 
is  a  sum  of  N,  equispaced  impulses: 

w, (t)  =  A.  T'  w,  6(t  -  mA.)  , 

1  lm  1  (24) 

where  {wlm}  is  a  finite  length,  real  sequence,  symmetric  about  m  =  0;  this  corres¬ 
ponds  to  discrete  sampling  of  waveform  x(t)  at  time  spacing  A,.  The  p-th  piece  of 
weighted  data  is,  from  (13), 

y  (t)  =  x(t)  At  £  wlm  6  ft  -  j-  -  pS  -  n.At) 
m  f 

•  4t  £  V  4  -  T-  -  Ps  -  "0  •  <25> 

til  / 


where  weighted  sample 
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The  first-stage  correlation  estimate,  (16),  becomes 


A 

Rj(t)  = 


E 

k 


Ik 


6(t  - 


kAt) 


where  the  area  of  each  impulse  is  given  by 


A 

R 


D  (p) 

Rk 


(27) 


(28) 


and 

=  A  V  V  V 

k  t  ^  ^p,m  ^p,m-k 

a  (29) 

The  last  quantity  is  the  sample  correlation  of  the  p-th  set  of  samples,  and  &lk  is  their 
average  over  the  total  of  P  pieces. 

The  first-stage  spectral  estimate  is  the  Fourier  transform  of  (27)  as  usual: 


A  A 

G1  Cf )  =  At  2-  Rik  exp(-i2irfkAt) 


(30) 


which  is  finite  for  all  f  and  is  of  period  1/A,  in  f.  An  alternative  expression  is 
available  by  substitution  of  (25)  in  (14): 


A 

Gj(f) 


1  P'1 

=  fE 

p=0 


£  e*p(-i2*fmV 

m  F 


(31) 


These  two  expressions  hold  for  arbitrary  f;  either  one  can  be  used  to  obtain  the 
first-stage  spectral  estimate.  If  we  restrict  our  calculations  of  interest  to  multiples  of 
some  frequency  increment  AF,  (31),  for  example,  specializes  to 

A  1  P'1  1 2 

Gi(qAF)  =  p  £  A  £  ypm  exp(-i2™qAtAF)  , 

P=0  m  '  (32) 


where  q  is  an  integer.  At  this  point,  there  needn’t  be  any  relation  between  A,  and  AF; 
we  can  calculate  the  spectral  estimate  at  any  frequencies  we  please.  However,  a 
favorite  choice  for  computational  purposes  is  to  choose  frequency  increment 


A 


F 


N  =  power  of  2  , 


to  get  the  special  digital  processing  result 


At  EyPm  exP(-i2irm9/N) 
m  * 


2 


(33) 


(34) 
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which  is  recognized  as  the  familiar  power  average  of  N-point  FFTs  of  weighted  data 
sequences.  All  the  impulsive  functions  in  (24)-(27)  have  dropped  out  of  first-stage 
spectral  estimates  (30)-(34). 

The  temporal  window  associated  with  impulsive  weighting  (24)  is  its  Fourier 
transform 


Wj(f)  -  A  £  w  exp(-i2irfmAt) 
m 


(35) 


Having  picked  an  impulsive  temporal  weighting  w,(t),  we  are  still  free  to  select  the 
lag  weighting  or  lag  window  as  we  please.  For  example,  for  any  lag  weighting 
function  w2(r),  (19)  and  (27)  yield  second-stage  correlation  estimate 


A 

R, 


(t)  =  A 


z 

k 


w2(kV  Rlk 


6(t  -  kAt)  =  At 


z 

k 


2k 


6  (t  *•  kAt) 


(36) 


The  corresponding  second-stage  spectral  estimate  is  the  Fourier  transform 

G2(f)  =  At  £  w2(kAt)  Rlk  exp(-i27tfkAt)  =  At  £  R2k  exp(-i2irfkAt) 
k  k 

(37) 

which  is  everywhere  finite  and  has  period  1/A,  in  f.  Evaluation  of  (37)  can  therefore 
be  confined  to  |f|  <  (2A,)"1. 

These  results  apply  for  general  lag  weighting.  A  specific  choice  is  the  lag  window 
with  N2  equispaced  nonzero  impulses: 


W2(f)  =  Af  £  W2n  6(f  -  nAf)  • 


(38) 


Frequency  spacing  Af  need  not  be  related  to  time  spacing  A,  in  (24),  nor  to  frequency 
increment  AF  used  in  the  frequency  and  FFT  calculations  above  in  (32)-(34).  Also 
there  are  no  relations  between  the  real  symmetric  sets  of  numbers  {wlm}  in  (24)  and 
{W2n}  in  (38).  Substitution  of  (38)  in  (18)  yields  for  the  second-stage  spectral 
estimate 

G2(f)  = Af  £  W2n  Gl^f  "  nAf^  » 

n  (39) 


which  is  a  local  average  (of  the  first-stage  estimates)  in  the  band  about  the 
frequency,  f,  of  interest.  Equation  (39)  is  a  discrete,  quadratic,  frequency¬ 
smoothing  operation.  In  fact,  (39)  holds  for  lag  window  (38)  and  any  temporal 
weighting  w,(t);  it  is  not  limited  to  the  discrete-time  form  (24). 

If  we  limit  our  calculations  of  <^(f)  to  multiples  of  frequency  increment  AF  as  in 
(32),  then  (39)  yields 

G2(qAF)  =  Af  £  W2n  G^qAp  -  nAf)  ;  (4Q) 
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we  can  use  (30)  for  on  the  right-hand  side.  Finally,  if  we  take  frequency  in¬ 
crement  (33)  and  frequency  spacing  Af  =  (MA,)'1,  where  integer  M  is  a  submultiple 
of  N,and  MA,  is  of  the  order  of  2L,,  the  FFT  results  of  (34)  can  be  employed  in  (40). 
More  will  be  said  later  on  the  choice  of  frequency  spacing  Af. 

The  variety  of  forms  available  at  different  stages  of  the  data  processing  illustrates 
a  great  deal  of  flexibility  in  exactly  how  the  available  data  x(t)  is  processed.  For 
example,  one  might  first  evaluate  via  FFT  procedure  (34).  Then,  since  (30)  can  be 
expressed  as 

Si(n^)  ■  \  £  “lk  .  (4|) 


it  follows  that  the  complete  nonzero  portion  of  correlation  sequence  {ftlk}  is 
recoverable  from  the  set  of  numbers  {6,(-£j-  )}^->  if  N^2N,-1,  where  N,  is  the 
number  of  nonzero  weights  {w,k}  in  (24)  (see  reference  18).  On  the  other  hand,  for 
N<2N,-1,  the  inverse  FFT  of  {^(-^j-)}^-1  would  yield  &)k  only  for  |k|  <  N-N, 
(reference  18,  eq.  (15));  thus  the  central  values  of  A  lk  are  recoverable  from  6,.  Then 
second-stage  correlation  estimate 


R2k  w2(kV  Rlk 


(42) 


follows  from  (36),  and  the  final  spectral  estimate  follows  from  (37).  The  lag 
weighting  samples  in  (42)  are  arbitrary;  thus  this  is  a  very  general  procedure  for 
obtaining  estimate  62(f)  at  any  f. 

The  relations  in  this  subsection  hold  for  arbitrary  values  of  t,  f,  and  q.  However, 
the  functions  of  t  are  impulsive,  and  are  zero  outside  limited  ranges,  while  the 
functions  of  f  and  q  are  periodic.  These  properties  should  be  utilized  in  any  com¬ 
puter  processing  technique  employing  these  forms.  Some  further  useful  properties 
and  interrelationships  of  the  sampled  lag  weightings  and  lag  windows  are  presented 
in  appendix  B. 
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Average  Value  of  Spectral  Estimate 

We  now  return  to  the  general  situation  for  both  the  temporal  and  the  lag 
weightings;  that  is,  we  do  not  presume  discrete  sampling  in  time  or  discrete 
smoothing  in  frequency.  From  (16)  and  (13),  the  mean  value  of  the  first-stage 
correlation  estimate  is 

Av{R1(t)  )  =  i  j  dt  x(t) x*  (t  -  t)  Wj  -  ~y  -  psj 

*  wl^  "  T  '  ~T  '  pS)  =  x^t^*(t  -  t)  ♦j(t)  =  R(t)  ^(t)  (43) 
where  R(t)  is  the  true  correlation  of  stationary  process  x(t),  and  where 

("O  =  f  dt  wt  (t)  w  (t  -  x) 

1  J  1  1  (44) 

will  be  called  the  correlation  of  real  temporal  weighting  w,(t);  see  the  footnote  to 
(9).  We  have  not  presumed  process  x(t)  Gaussian;  relation  (43)  holds  for  any 
stationary  process  x(t). 

Since  the  first-stage  spectral  estimate  6,(f)  is  a  linear  operation  (Fourier  trans¬ 
form)  of  ft  ,(t),  the  mean  value  of  6  ,(f)  is  the  Fourier  transform  of  (43);  that  is, 

AvfGjff)}  =J  dt  exp(-i27rfx)  R(t)  ^(t) 


=  G(f)  S  wf(f)  =  Jdu  Gj(u)  wj(f  -  u)  ,  (45) 

where  G(f)  is  the  true  spectrum  of  x(t),  i.e.,  Fourier  transform  of  R(t),  and  we  have 
Fourier  transformed  (44)  by  interchanging  integrals  and  using  temporal  window 


Wx(f)  =  J  dt  exp(-i2irft)  w^t) 


(46) 


The  convolution  result  in  (45)  is  a  familiar  one  for  the  standard  FFT  processing  of 
weighted,  overlapped  data  segments;  see  reference  9,  eq.  (5),  for  example.  Window 
W,(f)  is  real  and  even  about  f  =  0,  since  weighting  w,(t)  is  real  and  even  about  t  =0. 

The  mean  value  of  the  second-stage  correlation  estimate  follows  immediately 
from  (19)  and  (43): 

Av(R-(t)  }  =  R(t)  4>t  ( t)  w.(t)  =  R(t)  w  (t)  , 

*  1  L  “  (47) 


where 
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we(x)  =  ^(x)  w2(x)  (48) 

is  called  the  effective  (overall)  weighting  of  the  generalized  spectral  analysis 
technique.  It  incorporates  the  temporal  weighting  w,(t)  through  its  correlation  <^,(t) 
defined  in  (44),  and  it  involves  lag  weighting  w2(t)  directly.  Reference  to  (44)  and  to 
figure  2,  for  a  typical  temporal  weighting,  shows  that  $,(t)  is  as  depicted  in  figure  4; 
4|(t)  extends  over  (-L,,  L,)  and  is  zero  for  |t|  >  L,.  Since  the  effective  weighting 
we(r)  in  (48)  involves  $,(t)  and  w2(t),  we  now  see  the  reason  for  the  particular 
choices  of  L,  and  L2  in  figures  2-4.  Specifically,  ±L,  and  ±L2  measure  the  non¬ 
zero  extent,  in  the  r-domain,  of  the  functions  that  are  relevant  to  the  effective 
weighting.  Although  L,  measures  the  nonzero  extent  of  temporal  weighting  w,(t)  in 
the  time  domain  in  figure  2,  and  the  nonzero  extents  of  $|(t)  and  w2(t)  are  2L,  and 
2L2  in  figures  4  and  3,  respectively,  we  will  nevertheless  refer  to  L,  and  L2  as  the 
“lengths”  of  $,(t)  and  w2(t),  respectively,  in  the  t-domain,  for  convenience. 


Figure  4.  Correlation  $,(t)  of  Temporal  Weighting  w,(t) 


In  appendix  C,  ^(r)  is  evaluated  for  the  class  of  temporal  weightings* 

Wj(t)  =  ]£  ok  exp(i2irkt/Lj)  for  |t|  <  Lj/2  ,  ^ 

which  includes  a  wide  variety  of  weightings  such  as  rectangular,  Hamming, 
Hanning,  Blackman,  Harris,  and  the  recent  optimal  weightings  of  Nuttall,  reference 
19.  Specializations  to  real  symmetric  {ak}  and  to  a  limited  number  of  nonzero 
coefficients  are  also  made  in  appendix  C. 

Finally,  since  second-stage  spectral  estimate  <J2(f)  is  a  Fourier  transform  of  &2(t), 
its  mean  value  follows  from  (47)  and  (48)  as 


A 


Av{G2(f)}  =  G(f)  0  We(f) 


(50) 


where 


We(f)  =  fin  exp(-i2irfT)  v^Cx)  =  W^(f)  0  W2(f) 


(51) 


•For  brevity,  here  and  later,  we  omit  the  “0  otherwise”  statement  that  applies  for  |t|  >  L,/2,  as  was 
done  in  (3). 
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is  the  effective  (overall)  window  of  the  generalized  spectral  analysis  technique  of 
interest  here.  The  result  in  (51)  follows  by  Fourier  transformation  of  the  product  in 
(48)  and  use  of  (44)  (just  as  done  in  (45)).  Relation  (50)  is  a  simple  and  informative 
one  for  the  average  spectral  estimate;  it  enables  ready  determination  of  the  amount 
of  spreading  caused  by  particular  choices  of  temporal  and  lag  windows.  It  holds  for 
any  stationary  process  x(t)  with  spectrum  G(f);  thus  x(t)  needn’t  be  a  Gaussian 
process  for  (50)  to  hold  true. 

As  a  special  case  of  (50),  consider  lag  weighting  w,(x)  to  be  1  for  all  t.  Then 
W2(f)  =  <5(0  and  (50)  reduces  to  the  result  in  (45)  as  expected,  since  we  are  em¬ 
ploying  no  lag  weighting  at  all  in  this  case. 

As  another  special  case,  let  temporal  weighting  w,(t)  be  1  for  all  |t|  <  L,/2  and  let 
L,  =  T,  L,  «  T.  This  corresponds  to  Blackman-Tukey  processing.  Then  W^(f)  is 
proportional  to  sinc2(Tf),  which  is  much  narrower  in  f  than  W2(f),  meaning  that 
We(f) »  W2(f),  the  lag  window  alone. 

Interpretation  of  the  response  of  the  effective  window,  We(f),  via  convolution 
(51)  can  sometimes  be  deceiving,  and  it  may  be  helpful  and  necessary  to  resort  to 
(48).  For  example,  suppose  w,(x)  is  1  for  |t|  <  L2  and  0  otherwise,  where  L,  >  L,. 
Then  (51)  says  that  we  have  to  convolve  sinc(L2f),  which  has  -6.63  dB  side  lobes, 
with  Wj(f).  Our  first  impression  would  be  that  We(f)  is  bound  to  have  bad  side  lobes 
regardless  of  the  temporal  window.  But  recourse  to  (48)  and  figure  4  immediately 
reveals  that  we(x)  =  $((x)  for  all  x,  and  that  W2(f)  is  totally  irrelevant,  provided  that 
L2  >  L,.  The  scaling  of  $,(x)  by  a  constant  in  (48),  over  the  range  of  nonzero  $,(x), 
obviously  has  no  effect  on  the  relative  side  lobes  of  We(f).  Furthermore,  the  actual 
calculation  of  the  effective  window  via  (51)  is  often  tedious,  whereas  a  Fourier 
transformation  of  the  product  in  (48)  is  a  reasonable  approach,  even  if  only  by  an 
FFT. 
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Constraint  on  Temporal  and  Lag  Weighting  Lengths 

The  effective  window  We(f)  was  presented  in  (51).  Its  “width”  is  given  ap¬ 
proximately  by  the  sum  of  the  widths  of  the  temporal  and  lag  windows.  As 
discussed  earlier,  we  wish  to  constrain  the  effective  bandwidth  Be  of  We(f),  so  as  to 
be  able  to  fairly  compare  the  performance  of  different  spectral  analysis  techniques. 
The  effective  bandwidth  is  given  by  (1)  and  can  be  developed  as 

[Jdf  We(f)]2  w2(0)  *2(0)  w2(0) 

e  Jdf  W2(f)  Jdr  w2(t)  Jdr  *2(t)  w2(t) 


where  we  have  used  Parseval’s  theorem,  the  Fourier  transform  relationship  in  (51), 
and  (48).  Since  Be  is  to  be  considered  fixed,  (52)  forces  a  relationship  between 
lengths  L,  and  L2  of  $,(r)  and  w2(t). 

For  example,  consider  rectangular  temporal  and  lag  weightings  (this  is  not  a 
practical  case  and  is  presented  only  for  illustration  purposes): 

Wj(t)  =  1  for  | t |  <  Lj/2  , 

4>1  Ct)  =  Lj  -  |  t  |  for  |t|  <  Lj  , 


w2(t)  =1  for  I t I  <  L2 


Then  (52)  yields 


where 


Lm  =  rain(Lj ,  L2) 


(53) 


(54) 


(55) 


Given  a  value  of  Be,  (54)  can  be  considered  as  an  equation  for  L2  in  terms  of  L,,  or 
vice  versa.  Here  we  have  fixed  the  shapes  of  the  weightings  and  are  varying  the 
lengths  so  as  to  realize  the  specified  frequency  resolution  Be. 
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Generally,  the  qualitative  character  of  the  interrelationship  between  lengths  is  as 
depicted  in  figure  5,  for  fixed  weighting  shapes.  The  larger  one  of  the  lengths  is 
made,  the  smaller  the  other  length  can  get  and  still  satisfy  the  specified  bandwidth 
of  the  effective  window.  For  a  specified  pair  of  shapes  for  w,(t)  and  w2(t),  a  plot 
like  figure  5  can  be  used  in  two  different  ways.  If  we  pick  a  value  for  BeL,,  this 
determines  BeL2  and  hence  L2/L,.  On  the  other  hand,  choice  of  a  value  for  the  ratio 
of  lengths,  L2/L,,puts  a  line  through  the  origin  of  slope  L2/L,,  and  thereby 
determines  BeL,  and  BeL2  where  the  line  intersects  the  curve.  We  note  therefore  that 
knowledge  of  one  of  the  following  three  quantities  determines  the  other  two:  L2/L,, 
BeL„BeL2. 

The  limiting  parameter  values  on  figure  5  are  determined  as  follows:  as  L,  -*  00 , 
then  L2  -*  L2(min),  where  now  (from  (52)) 


e 


w2(t) 


-2 


W2(0) 


2L2(min)  c{*2)  . 


(56) 


Here,  c{  }  is  a  “shape  factor”  defined  for  any  limited-duration  function  g  as 

,2 

c{  '  i  i  ■  i eix) 

ZL 

g 


;{g,  -  jf-  /  d,  [|W]  . 


(57) 


where  it  is  assumed  that 


g(x)  =0  for  | r |  >  Lg  ; 


(58) 


Figure  5.  Interrelationship  of  Lengths  L(  and  L2  for  Fixed  Shapes 
of  the  Temporal  and  Lag  Weightings 
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that  is,  2Lg  is  the  nonzero  extent  of  g(-r).  Shape  factor  c{g}  is  independent  of  the 
magnitude  of  g(-r)  and  of  its  length  on  the  T-scale.  Thus  (56)  yields 
BeL,(min)  =  (2  c{w2})-*,  which  is  entered  on  the  ordinate  in  figure  5.  However, 
since  L,  is  limited  by  T,  the  dotted  portion  of  the  curve  on  the  far  right  is  not  at¬ 
tainable. 


Conversely,  if  instead,  L2  -*  00 ,  then  L,  -*■  L^min),  where  now  (from  (52)) 


(t) 


=  2L1  (min)  c { 4>  } 


(59) 


that  is,  BeL,(min)  =  (2  c{$,})_1.  This  value  is  attainable;  it  corresponds  to  no  lag 
weighting.  The  ratio,  L2/L,,  of  weighting  lengths  can  take  values  in  the  range 
L2(T)/T  to  for  BeT  »  1 ,  this  constitutes  the  range  from  almost  zero  to  infinity. 

Since  the  shape  factors  c{w2}  and  c{ <!>, }  are  important  limits  on  the  weighting 
lengths,  tables  of  their  numerical  values  for  a  number  of  useful  weight  functions  are 
given  below.  The  weightings  listed  under  C5,  C3,  Cl  are  those  given  in  reference  19, 
figures  10,  11,  12,  respectively;  the  notation  means 

C5:  continuous  fifth  derivative  of  weighting 
C3:  continuous  third  derivative  of  weighting 
Cl:  continuous  first  derivative  of  weighting 


For  the  class  of  lag  weightings  given  by 


w9(t)  »  £  a,  cosCirkr/L.)  for  |t|  <  L.  , 
k>.0  K  2  2 


(61) 


the  shape  factor  is 
c{v 

This  is  evaluated  for  several  weightings  in  table  1 . 


rdT  h(I)T .  “o  <?(ai  ♦  4  »  •••) 

2>  21,  )  [»2(0,J  (.  ..7  • 


(62) 


Table  1.  Shape  Factor  for  Lag  Weighting  w2(r) 


W2(t) 

c{w2(r)} 

Rectangular 

1.000 

Hanning 

.3750 

Hamming 

.3964 

Blackman 

.3046 

C5 

.2256 

C3 

.2442 

Cl 

.2558 

18 


TR  6459 


For  the  class  of  temporal  weightings  given  by 

w,(t)  =  £  av  cos(2nkt/L1)  for  |t|  <  L./2 
1  k>0  k  1  1  (63) 

the  correlation  4i(t)  is  evaluated  in  general  in  appendix  C.  The  shape  factor  of  $,(r), 


cUt)  = 


*j(T) 


(64) 


can  then  be  evaluated  numerically  and  is  given  in  table  2. 

Table  2.  Shape  Factor  for  Correlation 
$j(t)  of  Temporal  Weighting  w,(t) 


W,(t) 

c{<M 

Rectangular 

.3333  = 

1/3 

Hanning 

.2405  = 

(8n2  +  35)/(48n2) 

Hamming 

.2628 

Blackman 

.2073 

C5 

.1545 

C3 

.1678 

Cl 

.1763 

Plots  of  the  relationship  between  L,  and  L2  dictated  by  (52)  are  given  in  figure  6 
for  various  combinations  of  temporal  and  lag  weightings.  For  a  rectangular  lag 
weighting,  the  curve  will  actually  reach  BeL,  =  (2c{ 4, })-1  when  L2  =  L,;  then  the 
curve  goes  vertically  up  from  this  point  for  L2  >  L,  (see  figure  6A).  The  procedure 
for  the  evaluation  of  figure  6  is  as  described  under  figure  5;  namely,  pick  a  value  for 
L2/L,,  compute  BeL,  via  (52),  and  then  compute  BeL2  =  BeL,*(L2/L,). 

If  the  maximum  segment  length,  L,,  is  specified  (as  for  example,  when  the 
maximum  FFT  size  and  the  time-sampling  spacing  A,  are  fixed),  under  what  con¬ 
dition  can  a  desired  effective  frequency  resolution,  Bd,  be  met?  The  answer  to  this 
question  is  available  from  figure  5;  namely,  we  see  that 


BeLj  >  (2c{«p1» 


or 


B 

e 


~  2c{^1)L1 


(65) 


Thus  if  desired  resolution  Bd  is  greater  than  or  equal  to  the  right-most  term  of  (65), 
there  exists  a  choice  for  lag  length  L2  that  will  yield  the  desired  frequency  resolution. 
The  shape  factor  in  (65)  depends  only  on  the  temporal  weighting  w,(t). 
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6A.  Rectangular  Temporal  Weighting,  Rectangular  Lag  Weighting 


0  1  2  3  4  5  6 


6B.  Rectangular  Temporal  Weighting,  Hanning  Lag  Weighting 


Figure  6.  Allowed  Lengths  of  Various  Temporal  and  Lag  Weighting  Pairs 
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Figure  6.  Allowed  Lengths  of  Various  Temporal 
and  Lag  Weighting  Pairs  (Cont’d) 
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6E.  Rectangular  Temporal  Weighting,  Cl  Lag  Weighting 


6F.  Hanning  Temporal  Weighting,  Hanning  Lag  Weighting 


Figure  6.  Allowed  Lengths  of  Various  Temporal 
and  Lag  Weighting  Pairs  (Cont'd) 
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Normalization  of  Weightings 

The  average  of  the  first-stage  correlation  estimate  was  given  in  (43).  For  t  =  0,  it 
yields 


Av{Rj(0)}  = 


R(0)  *,(0) 


(66) 


Since  R(0)  is  the  true  power  in  the  process  x(t)  under  investigation,  it  is  convenient 
to  normalize  according  to 


1  =  ^(0)  =  f  dt  Wj(t) 


(67) 


Then  A,(0)  is  an  unbiased  estimator  of  R(0). 

Additionally,  from  (47)  and  (48),  we  have,  for  the  second-stage  correlation 
estimate, 


Av{R2(0)}  =  R(0)  4^(0)  w2 (0) 


(68) 


Therefore,  in  addition  to  (67),  we  also  set  lag  weighting  value 


w2(0)  =  1 


(69) 


making  ft2(0)  an  unbiased  estimator  of  R(0).  There  follows,  for  the  effective 
weighting, 


w  (0)  =  1 
e 


(70) 


Since  there  is  no  significant  loss  of  generality,  the  normalizations  in  (67),  (69),  and 
(70)  will  be  used  in  the  rest  of  this  report. 


Discrete-Time  Processing 

For  the  impulsive  temporal  weighting  introduced  in  (24)-(29),  the  normalization 
(67)  must  be  modified  somewhat,  since  the  integral  of  wf(t)  in  (24)  would  be  infinite. 
We  resort  to  (28)  and  require  that  the  origin  value  of  the  sample  correlation  satisfy 
the  unbiased  requirement  that 


Av{R10)  =  R(0) 


(71) 


Reference  to  (28),  (29),  (26),  and  (43)  yields 

AV(«10>  ■  A'',JoP)> '  E  V5  *  E  “L  M0) 

m  r  m 


(72) 
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therefore  the  normalization  is 


(73) 


This  is  the  discrete  analog  to  the  integral  constraint  in  (67). 

The  correlation  4,(t)  of  temporal  weighting  w,(t)  in  (24)  is  given  by  (44)  as  usual 
and  is  expressible  as 


where 


4>1  Ct)  =  At  <|>lk  6(t  -  kAt) 

T"*  w,  w  , 

*ik  t  l,m  l,m-k 

m 


(74) 


(75) 


Thus  we  see  that  (73)  is  tantamount  to  4io  =  1,  which  replaces  the  constraint  4,(0)  =  1 
in  (67)  for  the  continuous  temporal  weighting  case. 


If  we  also  require  that  second-stage  correlation  estimate  (42)  satisfy  the  unbiased 
requirement 


Av(R20>  =  R(0) 


(76) 


then,  as  before,  we  require 


w2(0)  =  1  . 


(77) 


Finally,  the  effective  weighting  becomes,  upon  use  of  (74), 

w  (t)  =  $  (t)  w  (t)  =  A  £  <f>-.  w  (kA  )  6(t  -  kA  )  . 

e  1  *  tk  ^t  1  (78) 

The  normalizations  adopted  above  make  the  area  of  the  impulse  at  r  =  0  equal  to 

A  410  w  (0)  =  A  . 

X  (79) 
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Examples  of  Effective  Windows 

We  consider  first  a  rectangular  temporal  weighting  w,(t),  for  which  the 
correlation  is  triangular. 


=  1  -  l T I /Lj  for  |  t  |  <  Lj  , 

and  the  class  of  lag  weightings  as  given  earlier  by  (61): 

w-(t)  =  23  a,  cos(Trkr/L  )  for  |r|  <  L. 

k>_0  1 

Then  constraint  (52)  yields 


(80) 


(81) 


Be  L2  =  12 


(  dx  [i,  ak  co!("“)]  ['  ‘  v] 


V-l 


for  L2  < 


(82) 

The  effective  weighting,  wc(t),  is  given  by  the  product  of  (80)  and  (81);  its  Fourier 
transform  is  the  effective  window 


2L 

V«  •  -7  fn  a. 

it  k>0 


1  .  !ll \  2nv  sin(2irv)  -  - 

L1  /  4v  -  k 


+  -  (-l)k  cos(2ttv)]  j  —  —  k 


where 


(4/ .  7) 


2  (  for  L2  -  L1  »  (83) 


v  =  BLu  ,  u  =  f/B  . 

6  2  e  (84) 

Although  (82)  and  (83)  could  be  extended  to  the  case  where  L2  >  L|(  that  range  is  not 
of  practical  interest,  as  will  become  apparent  later. 

The  numerical  procedure  for  evaluation  of  the  effective  window  is  to  first  select 
the  shape  of  the  lag  window  by  specifying  coefficients  {ak}.  Then  we  choose  a  value 
for  L2/L,  and  compute  BeL2  from  (82).  We  can  then  employ  (83)  and  (84)  to 
determine  We(f).  Four  examples  are  given  in  figure  7,  where  we  have  plotted 


Vf)  t 
dB  s  10 ysr  «  be 


(85) 
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7A.  Effective  Window  for  Rectangular  Temporal  Weighting 
and  Hanning  Lag  Weighting 


7B.  Effective  Window  for  Rectangular  Temporal  Weighting 
and  C5  Lag  Weighting 


Figure  7.  Examples  of  Effective  Windows  for  Rectangular  Temporal  Weighting 


TR  6459 


dB 


The  curve  corresponding  to  L2/L,  =  0  (i.e.,  L,  =  °°)  is  that  for  the  lag  window 
alone. 

The  overriding  impression  of  the  plots  in  figure  7  is  that  the  effective  window  has 
poor  side  lobe  behavior  and  decay  unless  L2/L,  is  chosen  very  small.  That  is,  the 
poor  side  lobe  behavior  of  temporal  window  Wf(f)  enters  the  convolution  (51)  for 
We(f),  and  is  difficult  to  suppress,  even  by  choice  of  good  lag  windows.  It  would  be 
desirable  to  realize  the  bottom-most  figures  in  each  of  these  plots,  since  these  latter 
curves  have  good  side  lobes  and  decay;  a  procedure  for  accomplishing  this  goal  is 
presented  in  the  next  section. 

The  situation  is  significantly  improved  when  the  temporal  weighting  is  tapered. 
An  example  for  Hanning  temporal  and  lag  weightings  is  given  in  figure  8.  The 
bottom-most  curve  has  an  eventual  18  dB/octave  decay  because  $,(t)  has  a 
discontinuous  fifth  derivative  at  t  =  0,  which  is  not  compensated  by  w2(t).  ($,(t) 
also  has  a  discontinuous  fifth  derivative  at  t  =  ±L,,  but  this  is  converted  to  a 
discontinuous  seventh  derivative  for  we(t)  by  means  of  w2(t)  when  L2  =  L(.) 


Figure  8.  Effective  Window  for  Hanning  Temporal  Weighting  and 
Hanning  Lag  Weighting 
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Discrete-Time  Processing 

The  temporal  weighting  w,(t)  is  given  by  (24),  and  its  correlation  $,(t)  by  (74), 
where  coefficients  {$lk}  are  given  by  (75).  The  Fourier  transform  of  (24)  leads  to 
temporal  window 


Vf)  *  \  £  wlm  exp  (-i2TTfmAt) 

m  (86) 

which  has  period  1/At  in  f,  and  is  real  and  even,  since  weight  sequence  {wlm}  is  real 
and  even. 

For  a  general  lag  weighting  w2(t),  the  effective  weighting  is  given  by  (78).  The 
effective  window  We(f)  is  given  by  (5 1 )  as  the  Fourier  transformation  of  (78): 

We(f)  =  At  £  *lk  w2(kAt)  exp(-i27rfkAt)  ,  (g?) 

which  also  has  period  1  /At  in  f,  and  is  real  and  even,  since  lag  weighting  w2(t)  is  real 
and  even.  This  result  holds  for  any  lag  weighting  w2(t)  and  is  a  very  useful  form  for 
computing  We(f)  for  any  value  of  frequency  f.  The  convolutional  form  of  (51)  is  not 
very  useful  for  computing  We(f)  for  general  W2(f). 

As  a  special  case,  we  can  evaluate  (87)  at  particular  frequencies 
fn  =  nAF  =  n/(NAt)‘1 ,  as  in  (33)  and  (34): 

*  \  £  *lk  “2(kV  '*P(-i2™k/N)  .  (W 

which  can  be  accomplished  as  an  N-point  FFT.  We  should  choose  N  large  in  order 
that  (88)  be  capable  of  tracing  the  fine  detail  of  We(f).  This  is  an  attractive  and 
efficient  way  to  evaluate  the  effective  window. 

A  Special  Lag  Window  for  Discrete-Time  Processing 

The  result  in  (87)  applies  for  discrete  time  sampling  and  arbitrary  lag  weighting. 
We  now  specialize  to  the  lag  window  given  in  (38): 

w,(f)  =  A,  £  W  6 (f  -  nAf)  , 

2  f  n  2n  f  (89) 

where  sequence  { W2n}  is  real  and  even,  and  frequency  spacing  Af  need  not  be  related 
to  time  spacing  A,  in  (24).  Then,  via  inverse  Fourier  transformation  (20),  the 
corresponding  lag  weighting  is 

w  (r)  =  A,  V  W  exp(i2TrnA.T)  , 

2  t  n  2  t  (90) 

and,  in  particular,  sample  values 

W-(kA  )  *  A.  V  W.  exp(i2nnk  A- A  )  , 

1  n  r  C  (91) 
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which  can  be  used  in  (87)  to  evaluate  We(f)  for  any  frequency  f  whatsoever. 

An  alternative  form  to  (87)  for  calculation  of  We(f)  for  this  special  lag  window 
(89)  is  afforded  by  substitution  of  (89)  in  the  convolutional  form  of  (51): 

V«  ■  sf  E  W2„  Wllf  -  "V  •  ,92, 

The  temporal  window  W,(f)  is  given  here  by  (86).  Equation  (92)  is  an  attractive 
form  when  the  number,  N2,  of  nonzero  coefficients  { W2n}  is  small  and  W ,(0  can  be 
evaluated  in  closed  form.  In  fact,  (92)  actually  holds  for  any  temporal  weighting 
w,(t);  it  is  not  limited  to  the  discrete-time  forms  (24)  and  (86).  Equations  (87)  and 
(92)  are  duals  in  the  sense  that  (87)  applies  to  any  w,(t)  and  an  impulsive  w,(t), 
whereas  (92)  applies  to  any  W^f)  and  an  impulsive  W2(f).  Either  equation  can  be 
evaluated  at  any  f  of  interest. 

Our  first  example  is  rectangular  temporal  weighting;  from  figure  2,  (24),  and 
(73),  wlm  =  Ly1/2,  where  L,  =  N,A,.  Then,  from  (86), 


w‘Cf)  =  L1 


sin(irL^f) 


Nj  sinfirLjf/Nj)  j 


Li  yv) 


(93) 


For  the  lag  window,  we  take  impulsive  form  (89);  then  (92)  and  (93)  give 


».<«  ■  h  4f  £  "2„Vli£  ■  "Lt4f)  ■ 

n  1  (94) 

Two  important  choices  yet  to  be  made  are  L,Af,  the  relative  frequency  spacing  used 
in  frequency  smoothing,  and  the  set  of  coefficients  {W2n}.  For  Hanning  frequency 
smoothing,  the  latter  is 

11/2  for  n  =  0 

1/4  for  n  =  +1 
— 

0  otherwise 

The  effective  windows  for  L,Af  =  1/2  and  1  are  given  in  figure  9  for  N,  =  32. 
Window  (94)  is  even  about  f  =  0  and  has  period  1/A,  in  f;  hence  only  the  region 
0,  (2A,)-'  is  plotted  in  f.  The  window  in  figure  9A  has  no  deep  notches  since  the 
frequency  displacement  (spacing)  Af  =  (2L|)_l  causes  the  notches  to  be  filled  in;  the 
window  for  Af  =  Ly1  in  figure  9B  reinforces  the  notches  and  has  a  significant 
shoulder  near  f  =  1/L,.  Both  windows  have  slow  decay  with  frequency  and  do  not 
have  significant  rejection,  even  near  Nyquist  frequency.  Closer  spacings  than 
(2L,)_1  do  not  improve  the  decay  or  rejection  capabilities;  wider  spacings  than  Ly1 
generate  humps  in  the  effective  window.  The  bad  features  of  rectangular  temporal 
weighting  are  not  undone  by  Hanning  frequency  smoothing;  see  also  figure  7A. 
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Figure  9.  Effective  Window  for  Discrete-Time  Rectangular  Temporal 
Weighting  and  Hanning  Frequency-Smoothing 
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The  second  example  is  rectangular  temporal  weighting  with  N,  =  32  and  rec 
tangular  frequency  smoothing  over  five  frequency  samples;  i.e., 

1/5  for  [  n  j  <_  2 
0  otherwise 

The  effective  windows  for  L,Af  =  1/2  and  1  are  given  in  figure  10.  The  main  lobe 
humps  in  figure  10B  are  caused  by  the  displacements  of  W|(f)  according  to  (92). 
Both  windows  again  have  poor  decay  and  poor  rejection;  however,  the  main  lobe  is 
more  box-like  in  shape  than  previously. 

The  third  example  is  identical  to  the  previous  one  except  that  N2  =  1 1 .  The  plots 
in  figure  1 1  reveal  that  the  main  lobe  is  quite  box-like,  but  the  decay  and  rejection 
are  no  better  than  previous  cases.  According  to  (92),  we  are  merely  taking  the  poor 
side  lobes  and  decay  of  W^(f)  and  moving  them  about,  but  not  improving  them  in 
any  way. 

The  last  example  in  this  subsection  is  Hanning  temporal  weighting  with  no 
frequency  smoothing  at  all.  The  effective  window  for  Nj  =  32  is  simply  W|(f)  and  is 
plotted  in  figure  12.  It  has  the  familiar  -31.5  dB  peak  side  lobe,  a  rapid  decay,  and 
significant  rejection  capability. 


32 


34 


Figure  11.  Effective  Window  for  Discrete-Time  Rectangular  Temporal  Weighting 
and  Rectangular  Frequency-Smoothing  with  N2  =11 
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Lag  Reshaping  for  Desired  Effective  Windows 

From  (48),  the  effective  weighting  of  the  generalized  spectral  analysis  technique  is 
given  by 


We(T)  =  ^(t)  W2(X) 


(97) 


where  $,(t)  is  the  correlation  of  temporal  weighting  w,(t);  see  (44).  Now  suppose 
that  for  a  given  temporal  weighting  w,(t),  with  associated  correlation  4>j(t),  we 
choose  lag  weighting 


w,(x) 

w2(t)  =  for  |t|  ~l2  -L1 


(98) 


where  wd(r)  is  a  desirable  weighting  with  wd(0)  =  1  (in  keeping  with  (67)  and  (70)) 
and 

w  (x)  =  0  for  |t)  >  L2 

(99) 


Notice  that  L2  >  Lt  is  disallowed  in  (98)  since  <j>,(T)  =  0  for  |t|  >  Lv  Then  sub¬ 
stitution  of  (98)  in  (97)  yields 


*cCO  =  Wd(x)  »  WeCf)  =  Wdtf) 


(100) 


That  is,  the  effective  weighting  and  window  are  equal  to  the  desired  behavior.  We 
have  “undone”  the  effects  of  bad  side  lobes  in  temporal  window  W  j(f)  by  reshaping 
according  to  lag  weighting  w2(x)  in  (98).  (The  effect  on  the  variance  of  the  second- 
stage  spectral  estimate  <52(f)  be  considerec*  later.) 

To  see  how  much  can  be  accomplished  by  this  approach,  some  attainable  ef¬ 
fective  windows  that  can  be  realized  via  lag  reshaping,  for  continuous  rectangular 
temporal  weighting,  are  given  in  figure  13  for  the  largest  possible  value  of  L2, 
namely,  L2  ==  Lv  Superposed  on  the  window  Wj(f)  for  rectangular  temporal 
weighting  are  the  effective  windows  for  four  candidate  lag  reshapings,  for  L2  =  L,. 
These  are  the  narrowest  possible  effective  windows  for  a  given  L(.  The  first  one  in 
figure  13A  corresponds  to  an  effective  Hanning  weighting.  The  peak  side  lobe  is 
only  reduced  from  -13.3  dB  to  -15.7  dB,  and  the  asymptotic  decay  is  improved  to  9 
dB/octave  from  6  dB/octave.  The  main  lobe  width  is  only  slightly  broadened. 

Much  greater  improvements  in  side  lobe  behavior  are  possible  with  other  lag 
weightings,  and  are  illustrated  in  parts  {B)-(D)  of  figure  13.  They  illustrate, 
respectively,  peak  side  lobe  levels  and  decays  of:  -30.5  dB,  21  dB/octave;  -41.3  dB, 
15  dB/octave;  and  -46.7  dB,  9  dB/octave.  The  deeper  peak  side  lobe  is  realized  at 
the  expense  of  a  slower  asymptotic  decay.  They  all  have  about  the  same  main  lobe 
width.  The  C5,  C3,  Cl  weightings  were  introduced  and  explained  in  (60). 
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13A.  Reshaping  to  Achieve  Hanning  Weighting 


13B.  Reshaping  to  Achieve  C5  Weighting 


Figure  13.  Some  Attainable  Effective  Windows  via  Lag  Reshaping 
for  Rectangular  Temporal  Weighting  and  L2  =  L. 
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Figure  13.  Some  Attainable  Effective  Windows  via  Lag  Reshaping  for 
Rectangular  Temporal  Weighting  and  L2  =  L,(Cont'd) 
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Figure  13  illustrates  how  advantageous  the  reshaping  technique  can  be  in  terms  of 
peak  side  lobe  and  asymptotic  decay,  although  the  main  lobe  width  is  significantly 
increased.  In  fact,  the  peak  side  lobe  at  f  £  1.5/L(  for  the  rectangular  window  is 
really  not  suppressed, so  much  as  it  is  smeared  out;  however,  the  other  peaks  of 
W|(f)  for  | f)  >  2/L,  are  indeed  significantly  reduced.  Thus  reduction  of  leakage  via 
lag  reshaping  is  a  very  effective  method,  provided  that  we  accept  the  nearest  side 
lobe  of  the  temporal  window;  this  conclusion  is  in  contrast  to  reference  20,  page  57. 
These  general  conclusions  on  lag  reshaping  hold  also  for  temporal  weightings  other 
than  rectangular,  although  the  exact  degree  of  improvement  will  be  different. 


If  L2  is  chosen  less  than  L,,  the  effective  windows  in  figure  13  are  simply 
broadened  according  to  the  ratio  L,/L2.  The  peak  side  lobe  levels  and  asymptotic 
decay  are  unchanged,  but  the  main  lobe  width  is  increased.  Here  we  are  presuming 
L,  fixed  and  decreasing  L2. 

If  we  insist  that  the  combination  of  temporal  weighting  w,(t)  and  lag  reshaping 
w2(t)  in  (98)  have  effective  bandwidth  Be,  then  use  of  (52),  (99)-(100),  and  (57)-(58) 
yields 


B  = 
e 


wj(0) 


“d<°> 


f di  w^(x)  J< 


dx  wd(x) 


2L^-cTw^T 


where  c{wd}  is  the  shape  factor  of  wd(t)  (see  table  1).  Thus 


BeL2  sTcTw^T  f°r  L2^L1  ; 


(101) 


(102A) 


the  limits  on  L,  in  (102A)  follow  from  (98)  and  (15). 

A  plot  of  the  interrelationship  between  L,  and  L2  (introduced  in  figure  5)  is  shown 
in  figure  14  for  the  case  of  lag  reshaping.  The  reason  that  the  plot  is  flat,  in  contrast 
to  figures  5  and  6,  is  that  the  shape  of  w2(r)  now  changes  as  L,  changes.  This 
behavior  is  discernible  from  (98),  since  the  denominator  varies  while  the  numerator 
remains  fixed  according  to  the  selection  of  wd(-r)  and  its  associated  bandwidth- 
length  factor  (102A). 


If  the  maximum  segment  length,  L,,  is  specified  (as  for  example  when  the 
maximum  FFT  size  and  the  time-sampling  spacing  A,  are  fixed),  the  condition  under 
which  a  desired  effective  frequency  resolution,  Bd,  can  be  met  is  given  by  figure  14. 
Namely,  we  see  that 


B 


e 


1 

-  2  c{wd)  Lj 


(102B) 


Thus  if  desired  resolution  Bd  is  greater  than  or  equal  to  the  right-side  of  (102B), 
there  exists  a  choice  for  segment  length  L,  that  will  yield  the  desired  frequency 
resolution.  The  shape  factor  in  (102B)  depends  only  on  the  desired  weighting  wd(-r). 
(See  (65)  and  the  accompanying  discussion  for  the  case  where  lag  reshaping  is  not 
employed.) 
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Variance  of  Spectral  Estimate 

Up  to  this  point,  we  have  presumed  nothing  about  the  process  x(t)  except  that  it 
be  (second-order)  stationary  during  the  observation  interval  T.  Now  we  make  a 
couple  of  assumptions  about  the  process  in  order  to  obtain  manageable  expressions 
for  the  variance  of  the  second-stage  spectral  estimate  62(f),  at  frequency  f.  Our  first 
assumption  is  that  the  true  spectrum  G  of  x(t)  varies  slowly  in  the  neighborhood  of 
the  frequency  of  interest,  f.  More  precisely,  from  (50)  and  (70),  we  obtain,  for  the 
mean  spectral  estimate, 

Av{G2(f)}  =  J du  G(f  -  u)  Wg(u)  -  G(f)  f  du  We(u)  =  G(f)  , 

(103) 

where  we  assume  that  spectrum  G  is  relatively  constant  in  the  frequency  band 
(f-Be/2,  f  +  Be/2);  i.e.,  the  only  region  where  effective  window  We  in  (103)  is 
substantially  nonzero  is  in  the  range  (-Be/2,  Be/2). 

Our  second  assumption  is  that  x(t)  is  a  complex  Gaussian  process.  The  variance  of 
62(f)  is  developed  under  this  assumption  in  appendix  D,  culminating  in  the  exact 
result  in  (D-13): 

Var{G2(f)}  =  Jf  da  dg  G(a)  G(g)  |y(f  -  a,  f  -  g)|2  Qp(S(a  -  8))  , 

(104) 


where  window  convolution  function 

y(x,  y)  =  J  du  W2(u)  Wj(x  -  u) 

and  periodic  function 

n  fnt  -  rsinc^  V 
QP(U)  =  |_P 

The  variance  result  in  (104)  does  not  require  that  spectrum  G  vary  slowly  in  the 
neighbor  of  f;  the  result  utilizes  only  the  Gaussian  assumption  on  the  process  x(t). 
The  temporal  and  lag  windows  contribute  through  the  window  convolution  func¬ 
tion  y,  while  the  shift  S  and  number  of  pieces  P  appear  through  the  periodic  func¬ 
tion  Qp. 

When  the  assumption  regarding  a  slowly  varying  spectrum  G  in  the  neighborhood 
of  frequency  f  of  interest  is  also  invoked,  (104)  simplifies  to  forms  given  in  (D-20) 
and  (D-24);  the  latter  is  a  “weighting  domain”  version  of  the  variance: 

Var{G2(f)}  =  G2(f)  ±  £  fl  -  J  dr  w2(x)  *3(x,  pS)  , 

P"  "  (107) 


W*(y  -  u) 


(105) 


(106) 
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where 

=  fdt  Wj  (t  *  w^t  .  ^V-1)  "i(‘  *  H-1)  W1  (C  * 

"  V-T*  -u)  (108) 

is  a  third-order  correlation  of  temporal  weighting  w,(t).  The  form  (107)  is  very 
useful  if  $3  can  be  evaluated  in  closed  form.  An  “ambiguity  domain”  version  of  the 
variance  is  given  by  (D-20). 
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Quality  Ratio 


The  quality  ratio  for  spectral  estimation  was  defined  in  (4).  With  the  aid  of  (103) 
and  (107),  it  is  given  by 

Var{S,(f)} 

Q  =  “IT - 

Av  (G2(f)} 

=  p  pj.pl1  "  "^)/dT  W2(t')  *3(t’  pS)  '  (109) 

Since  the  smallest  possible  value  of  Q  is  (TBe)->  (see  (11)  et  seq.),  the  normalized 
quality  ratio  is 


TB  P-1  /  ,  ,  \  _ 

q  •  TBe  =  ~r  5  P  v1  ’  p  )  J dT  W2(t)  *3(t»  pS) 


NQR 


(110) 


This  quantity  can  never  be  smaller  than  unity. 


If  we  employ  (52)  and  the  normalizations  (67)  and  (69),  the  convenient  form 


NQR  =  T 


1_ 

P 


J dt  w2(t)  4>j( 


2 

dr  w2(t)  ♦jft.pS) 
t) 


(Ill) 


for  the  normalized  quality  ratio  is  obtained.  We  are  interested  in  the  behavior  of  the 
normalized  quality  ratio  for  different  choices  of  P,  S,  wt(t),  and  w2(t).  The  con¬ 
straint  of  a  fixed  effective  bandwidth  Be  has  been  injected  into  the  normalized 
quality  ratio  via  the  use  of  (52)  in  (111).  The  quantities  and  needed  in  (1 1 1)  are 
given  by  (44)  and  (108)  respectively. 


Before  we  embark  on  particular  cases,  some  general  observations  on  overlap 
(shift  S)  are  in  order.  For  a  minimum  normalized  quality  ratio  (minimum  variance) 
with  each  temporal  weighting  w,(t),  we  should  use  approximately  the  optimum 
overlap  as  derived  in  reference  9.  For  example,  Hanning  temporal  weighting  should 
be  employed  with  approximately  62  percent  overlap,  although  there  is  only  an  8 
percent  loss  in  stability  if  50  percent  overlap  is  used  for  convenience  (reference  9, 
tables  5  and  6).  There  is  no  point  in  considering  excessive  or  inadequate  overlap, 
since  this  leads  to  excessive  computational  effort  or  more  variance,  respectively. 
Inadequate  temporal  overlap  cannot  be  made  up,  in  terms  of  variance  reduction,  by 
any  amount  of  quadratic  smoothing.  This  can  be  seen  by  observing  that  poor  first- 
stage  correlation  estimates  ft,(T)  are  merely  multiplied  by  lag  weighting  w2(t),  and 
are  not  improved  statistically  in  any  way  for  |t|  <  L2;  those  estimates  for  |t|  >  L2  are 
discarded  by  the  lag  weighting. 


Some  related  work  on  the  effects  of  windowing  on  stability  is  given  in  references 
21  and  22.  However,  the  present  report  is  more  thorough  and  detailed  in  its 
treatment  of  the  problem  and  the  inclusion  of  a  bandwidth  constraint. 
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Special  Cases  of  Generalized  Spectral  Analysis  Technique 

This  section  will  consider  several  special  cases  of  the  normalized  quality  ratio  and 
show  how  some  earlier  results  are  obtained  as  limiting  cases.  The  next  section  will 
treat  the  generalized  spectral  analysis  results. 


One  Piece,  P  =  1 

When  only  one  piece  is  used  in  the  first-stage  spectral  estimate  6((f)  in  (14),  we 
have  a  generalized  version  of  the  Blackman-Tukey  approach,  in  that  the  data  x(t) 
are  weighted  by  w,(t)  prior  to  computing  the  sample  correlation;  see  (16)  and  (13). 
Also,  we  allow  length  L,  <  T  (although  we  soon  show  that  the  best  L,  is  equal  to  T, 
the  available  record  length).  From  (111), 

fdi  w? (x)  4>2(t) 

NQR1  =  T  4 - 2 - 2 - 

]  dr  w2(t)  (t)  (H2) 

where  (using  (108)) 

<t>2 s  ^3^T» 0)  =  fdt  wj(t  +  wj(t  -  j)  =  /dt  Wj(t)  w^(t  -  t) 

(113) 


is  the  correlation  of  the  squared  temporal  weight  function  wf(t). 


Now  if  L,  »  L2,  w2  is  much  narrower  than  or  $2.  In  that  case,  the  exact  shape 
of  w2  is  irrelevant,  and  (1 12),  (113),  (44),  figure  2,  and  Schwarz’s  inequality  yield 

fdt  w*(t)  T 


NQR 


1 


♦2(0) 

♦?(0) 


*  T 


Id 


[/dt 


2  — 


>  z —  for  L,  >>  L 


1 


(114) 


Equality  in  (1 14)  results  if  and  only  if  wf(t)  is  constant  for  |t|  <  L,/2;  furthermore, 
the  best  value  for  L,  is  then  its  largest  allowed  value  T  (see  (15)),  in  which  case  we 
have  Blackman-Tukey  processing  and 


NQR j (rectangular  w^)  m  1  for  T  =  Lj  >>  ly 


(115) 


This  result  agrees  with  reference  2,  section  B.8.  It  should  be  noted  that  L,  »  L2 
implies  BeL,  »  BeL2  >  .5/c{w2}  ~  1,  according  to  figure  5  and  table  1;  thus  stable 
estimates  result  in  this  case. 


Instead  of  rectangular  temporal  weighting,  consider  Hanning  weighting: 


(116) 
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Then 

♦  j(0)  =  /dt  w2(t)  =  1  ,  4>2 CO)  =  / dt  wj(t)  =  — - 

and  (114)  yields 


(117) 


NQRj^ (Hanning  w^) 


55  T 
18  Lx 


for  >> 


(118) 


The  best  L,  is  again  T;  however,  the  minimum  value  of  the  normalized  quality  ratio 
is  then  35/18,  which  is  twice  the  value  in  (1 15)  for  rectangular  temporal  weighting. 
This  is  due  to  the  squandering  of  the  edges  of  the  available  data  record  by  the  small 
values  there  of  Hanning  temporal  weighting. 

Now  instead  of  assuming  L,  »  L2,  let  us  reconsider,  for  general  w,,  Lp  L2,  the 
normalized  quality  ratio  (1 12).  Since  w,(t)  is  zero  for  |t|  >  L,/2,  we  have  from  (44), 


b  f  T ) 

(t)  =  f  dt  w  (t)  wx(t  -  t)  for  |x|  <  Lj 
a(t)  1 

where 

.  a(r)  =  max(-L^/2,  -  L^/2  +  x)  I 

j  for  | t |  < 

b(x)  =  min(Lj/2,  Lj/2  +  x)  j 

Then  by  Schwarz’s  inequality,  (113),  and  (120), 

j  rb(x)  j  j  /•  b(x) 

4j(x)  <J  dt  wj(t)  w*(t  -  xOj  dt  1 

a(x)  -a(x) 

=  4>2(x)  (Lx  -  |  x  | )  for  | t |  <  Lx  . 


(119) 


(120) 


(121) 


Equality  is  realized  in  (121)  if  and  only  if  w,(t)  is  constant  for  |t|  <  L,/2;  that  is,  the 
best  temporal  weighting  for  maximum  stability  is  rectangular  when  P  =  l.  This 
conclusion  holds  regardless  of  the  form  of  lag  weighting  w,(t)  or  the  relative  sizes  of 
L,  and  L2. 

As  an  example,  for  rectangular  temporal  weighting, 

Wj(t)  =  Lj”**  for  1 1 1  <  Lj/2  , 

♦  j(t)  =  1  -  M/Lj  for  | x  |  <  L,  , 

<t>2(x)  *  -f-C1  -  M/Li)  for  I T I  <  l2  ’  (,22) 
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and  (112)  yields  for  general  w2(t), 


T  C dT  (*  •  h-) 

NQR  (rect.  w  )  =  —  — - - 

1 


The  ratio  of  integrals  is  obviously  greater  than  1.  For  a  monotonically  decreasing 
lag  weighting  w2(x)  of  fixed  shape,  the  ratio  of  integrals  is  minimized  by  choosing-L, 
as  large  as  possible.  Since  the  leading  factor  also  has  the  same  behavior,  the  best 
value  for  the  normalized  quality  ratio  is 


NQRjCrect.  wJ)  = 


^  /T  dT  “2W(1  ' 

Jj dT  w2(t)(‘  -  ¥) 


for  L.  =  T 


We  cannot  give  numerical  values  to  this  ratio  of  integrals  until  we  select  a  lag 
weighting  w2(t)  and  determine  the  specific  value  of  L2(T);  see  figures  5  and  6.  But  if 
TBe  »  1,  which  is  the  usual  case  for  reasonably  good  spectral  estimates,  then 
L2«  L,  =  T  and 

NQRjfrect.  -  1  for  Lj  =  T,  TBg  »  1  . 

This  result  holds  independently  of  the  exact  shape  of  the  lag  weighting  w2(t);  thus 
we  could  choose  w2(t)  such  that  the  effective  weighting  we(r)  in  (48)  has  good  side 
lobe  behavior,  as  discussed  in  an  earlier  section. 

No  Quadratic  Frequency-Smoothing 

No  quadratic  smoothing  corresponds  to 

W2(f)  =  6(f), 

w2(x)  =  1  for  all  x  .  (126) 

Thus  L2  =  <»,  and  (109)  becomes 

Q  -  f  Ep  (i  -  J^L)  ♦st'.P3) 


(*  -  ¥)  • 
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since 


Jdz  <J>3(t,u)  =  J dr  j du  w1(u)w1(u  -  pJWjfu  +  t)  WjCu  +  t  -  y) 
=  J  du  Wjfu)  WjCu  -  p)  ^(y)  =  (J>j  Cw) 


(128) 


by  use  of  (108).  The  result  in  (127)  is  identical  to  reference  9,  equation  8,  when  we 
recall  definition  (109)  and  normalization  (67). 


Non-Overlapping  Segments 

Let  us  choose  time  shift  S  in  (13)  equal  to  the  segment  length  L, ;  this  leads  to 
abutting  time  segments.  From  (15),  we  have 


PLj  =  T 


(129) 


where  we  have  chosen  to  use  up  all  of  the  available  data  length.  (This  is  different 
from  the  earlier  subsection  for  P  =  1  where  we  allowed  L,  <  T.)  The  general 
normalized  quality  ratio  in  (1 1 1)  reduces  to 


fdi  w^(t)  <f>2  Cx) 

NQRp  =  L1  ~f - 2 - 2 -  * 

P  1  J  dx  (t)  *j(t) 


(130) 


where  we  used  (129),  the  fact  that  w,(t)  is  of  length  L,,  and  (113). 

Once  again,  we  refer  to  bound  (121)  and  the  fact  that  equality  is  realized  only  for 
a  flat  weighting  w,(t).  Thus,  from  (130)  and  (122), 


NQRp(rect.  w^) 


(131) 


for  any  (real  symmetric)  lag  weighting  w2(t).  The  ratio  of  integrals  is  obviously 
always  greater  than  unity;  therefore,  for  a  monotonically  decreasing  lag  weighting 
w2(t)  of  fixed  shape,  values  of  L,  large  in  comparison  with  L2  are  preferred. 
However,  L,  »  L2  means  that 


TB 

e 


PL,  B  >>  PL.  B  > 
1  e  2  e  — 


P 

2c{w2> 


P 


(132) 


according  to  (129)  and  figure  5.  Thus  large  time-bandwidth  products,  TBe,  are 
required;  also  P  must  be  kept  small  enough  to  realize  L|  »  L2.  In  this  case,  we  have 


NQRp(rect.  w^)  «  1  for  TBe  >>1,  Lj^  >> 


(133) 
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regardless  of  lag  weighting  w2(t).  Qualitatively,  when  L,  »  L2,  the  edge  effects  of 
segmenting  the  data  x(t)  are  negligible,  since  only  a  small  fraction  of  the  utilized 
correlation  values  that  can  be  calculated  from  a  record  of  length  T  are  neglected 
when  using  segments  of  length  L,.  Stated  alternatively,  all  the  first-stage  correlation 
estimates  that  are  used,  namely  R,(x)  for  |t|  <  L2,  have  the  same  quality  (stability) 
when  L2«  Lr 

An  example  of  the  exact  normalized  quality  ratio  for  Hanning  lag  weighting  and 
rectangular  temporal  weighting  is  afforded  by  substituting  the  equation 


w2(t)  =  cos2  for  M  <  L2 

'  2  > 


(134) 


into  (131)  (see  figure  7A  for  the  effective  window): 


NQRp(rect .  Wj,  Hann.  = 


fy '  4  /  2  L2 

I  dx  cos  x  II  -  —  j —  x 

J0  \  L1 


/. 


.  4/2  L2  f 

0  dx  cos  x  v  -  7177 


(135) 


where 


y  =  f  min/1,  . 


(136) 


Equation  (135)  is  plotted*  as  the  top  curve  in  figure  15 A.  As  expected,  the 
normalized  quality  ratio  tends  to  1  as  L2/Lj  tends  to  zero.  But  even  for  as  large  a 
value  as  L2/Lj  =  .5,  the  normalized  quality  ratio  has  increased  only  by  12  percent. 
Thus  the  penalty  in  increased  variance, for  not  realizing  a  small  ratio  for  L2/L„is  not 
severe. 

Also  plotted  in  figure  15A  is  the  normalized  quality  ratio  for  the  three  lag 
weightings  introduced  in  (60)  et  seq.  They  all  lead  to  smaller  values  of  the 
normalized  quality  ratio,  for  the  same  value  of  L2/Lp  in  fact,  lag  weighting  C5 
incurs  only  a  7  percent  increase  in  variance  when  L2/Lj  =  .5,  in  relation  to  the  ideal 
value  1.  The  reason  that  the  normalized  quality  ratio  is  lower  is  due  to  the  fact  that 
the  lag  weightings  drop  to  zero  faster  within  their  length  L2. 


Non-Overlapping  Segments;  Lag  Reshaping 


The  possibilities  of  lag  reshaping  have  been  discussed  earlier  with  regard  to  the 
mean  of  the  spectral  estimate  and  the  effective  window.  We  now  want  to  see  what 
effect  lag  reshaping  has  on  the  normalized  quality  ratio  in  (130).  Substitution  of  (98) 


in  (130)  yields 


NQRp(lag  reshaping) 


/ 


dx 


w2(x) 


4>2(t) 

♦  l(T) 


/ 


dx 


W?(T) 


for 


L2  ~ 


037) 


•The  quantity  TBe  is  not  involved  in  figure  15;  some  related  computational  considerations  are 
discussed  in  appendix  E. 
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Figure  15.  Normalized  Quality-Ratio  for  Rectangular  Temporal 
Weighting  and  No  Overlap 
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The  special  case  of  rectangular  temporal  weighting  is  obtained  by  employing  (122) 
in  (137): 

{  *  “>>  (>  - 

NQR  (rect.  w..  ,  lag  reshaping)  =  - ; -  for  L_  <.  L. 

p  i  l2  i  i 

{  dT  V°  038) 

The  division  by  ^(t)  in  reshaping  (98)  increases  the  variance  (for  a  specified  L2/Lj 
and  for  wd(T)  =  w2(t))  above  that  in  (131),  because  we  are  more  heavily  weighting 
regions  where  the  denominator  in  (138)  is  smaller. 

Equation  (138)  is  plotted  in  figure  15B  for  desired  effective  weightings  of  Han¬ 
ning,  Cl,  C3,  and  C5.  Notice  that  the  abscissa  is  now  limited  to  L2/Lj  <  1.  As 
expected,  the  normalized  quality  ratio  tends  to  1  as  L2/Lj  tends  to  zero;  that  is,  we 
can  do  lag  reshaping  for  good  side  lobe  behavior  and  lose  little  in  terms  of  stability, 
provided  that  Lx  is  chosen  sufficiently  larger  than  L2.  Of  course,  the  normalized 
quality  ratio  values  in  figure  15B  are  larger  than  those  in  figure  15A,  for  the  same 
value  of  L2/Lj.  As  an  example,  for  desired  effective  weighting  C,,  if  we  take 
=  2L2,  the  increase  in  variance  over  the  ideal  value  is  only  9  percent.  Thus  lag 
reshaping  is  an  attractive  procedure  for  spectral  estimation;  recall  from  figure  14 
that  L2  is  set  by  the  specified  Be  and  the  shape  of  wd(T). 
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General  Results  on  Stability 

We  now  return  to  the  general  normalized  quality  ratio  in  (111)  and  recall  con¬ 
straint  (15).  We  will  select  time  shift  S  according  to 

s  =  qLl  .  (139) 

where  q  =  q{w,}  is  a  fraction  specified  to  be  in  the  range  (0,  1]  and  is  dependent  on 
the  particular  temporal  weighting  w,(t)  employed.  The  observations  made  in  the 
paragraph  following  (111)  are  relevant  in  this  regard.  For  example,  with  no 
quadratic  frequency  smoothing  and  with  Hanning  temporal  weighting,  q  =  .39  (61 
percent  overlap)  is  virtually  optimum,  although  q  =  1/2  loses  only  8  percent  in 
variance-reduction  capability  (reference  9,  tables  4-7) .  We  also  select  equality  in  ( 1 5) 
so  as  to  make  maximum  use  of  the  available  record  length,  i.e.,  minimum  variance 
of  the  spectral  estimate.  Then  we  have 


(P  >.  1) 


(140) 


Thus  for  a  given  T  and  specified  shift  fraction  q,  L,  can  take  on  only  a  discrete  set 
of  values. 

Arbitrarily  large  values  of  P  are  not  allowed  in  (140),  because  this  would  result  in 
such  small  values  of  L,  that  the  bandwidth  constraint,  Be,  could  not  be  met.  From 
figure  5,  the  lower  bound  on  BeL,  limits 


2c  ( <t» « )  B  T  -  1 

p  =  1  +  - - - - - 

max  q 


(141) 


(Actually,  Pmax  must  be  the  integer  part  of  the  right-hand  side.)  Thus  Pmax  depends 
on  the  temporal  weighting  w,(t)  directly  through  its  shape  factor  c{i,}  and  in¬ 
directly  via  the  selected  shift  fraction  q  =  q{w,}.  For  q  =  1,  no  overlap,  (141) 
reduces  to  (E-l). 

When  P  =  Pmax,  L,  is  at  its  minimum  value,  and  L2  must  be  greater  than  T;  it  can 
be  <»  (see  figures  5  and  6A).  In  this  case,  there  is  no  quadratic  frequency  smoothing, 
and  we  have  the  situation  studied  in  detail  in  references  9  and  12,  and  mentioned 
earlier  in  ( 1 26)-(l 28).  At  the  other  extreme,  when  P  =  1,  we  have  Blackman-Tukey 
processing  with  the  generalization  that  the  temporal  weighting  need  not  be 
rectangular;  this  case  was  considered  in  the  previous  section.  The  range  of  values  of 
L2/L,  is  shown  in  (E-2)  and  (E-3)  to  be  very  wide  when  BeT  »  1,  which  is  a  usual 
practical  case. 

More  generally,  for  P  in  the  range  (1,  Pmax],  we  can  investigate  the  tradeoff 
between  the  amounts  of  temporal-  and  lag-weighting,  for  specified  resolution  Be 
and  for  specified  weighting  shapes  of  interest.  Below,  we  consider  the  two  cases  of 
rectangular  temporal  weighting  and  Hanning  temporal  weighting. 
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Rectangular  Temporal  Weighting 

It  was  shown  earlier  in  figure  15  that  rectangular  temporal  weighting  with  no 
overlap  results  in  small  values  for  the  normalized  quality  ratio,  whether  the  lag 
weighting  is  reshaped  or  not,  provided  that  L2  is  chosen  somewhat  smaller  than  L,. 
Now  the  question  arises  as  to  whether  one  should  use  any  overlap,  such  as  50 
percent,  with  rectangular  temporal  weighting. 

We  presume  BeT  »  1.  For  no  overlap,  the  estimate  of  first-stage  correlation 
R,(t)  at  t  ~  L,/2  has  only  half  the  degrees  of  freedom  as  the  estimate  at  r  =  0.  But 
with  50  percent  overlap,  the  degrees  of  freedom  for  estimation  at  t  =  L,/2  are  about 
the  same  as  at  t  =  0.  This  is  why  50  percent  overlap  for  rectangular  temporal 
weighting  appears  attractive. 

However,  for  estimation  of  R,(t)  at  t  >  L,/2,  we  still  do  not  get  as  many  degrees 
of  freedom  as  for  Blackman-Tukey  processing,  because  some  data  points  never 
interact.  For  example,  although  at  t  =  3L,/4  we  have  doubled  the  degrees  of 
freedom  by  using  50  percent  overlap,  we  still  have  only  about  half  of  the  number 
that  are  available  at  this  t  value  via  Blackman-Tukey  processing. 

In  order  to  ascertain  quantitatively  the  merit  of  overlapping  for  rectangular 
temporal  weighting,  we  have  evaluated  the  normalized  quality  ratio  (111)  for  lag 
reshaping  to  realize  a  desired  effective  weighting  equal  to  Cl  as  given  in  (60)  and 
(61)  (reference  19,  figure  12).  That  is,  in  (111),  we  use 

wd(T) 

w2Ct)=?7(TT  for  I T I  <  L2  ~  L1  * 

1  (142) 

where  $|(t)  is  given  by  (122).  In  addition,  we  need  the  third-order  correlation  (108), 
which  is 

$3(t,u)  =  -jj-  -  illj— liil\  for  | t  |  +  |u|  <  Lj 

1  1  (143) 

for  rectangular  temporal  weighting. 

For  the  two  cases  of  BCT  =  100andBeT  »  1000,  the  normalized  quality  ratio  has 
been  evaluated  for  q  =  1,  .75,  and  .5,  and  plotted  in  figure  16.  The  explanation  of 
the  behavior  of  the  curves  is  as  follows: 

q  =  1  No  Overlap 

If  P  =  I,  then  L2 «  L,  and  it  follows  that  for  all  |r|  <  L2,  R,(t)  is  estimated  with 
virtually  the  same  degrees  of  freedom  as  at  t  =  0,  where  we  have  the  maximum 
degrees  of  freedom  possible  to  estimate  R,(0).  As  P  increases  toward  Pmax,  then  L, 
tends  to  L2.  Now  R,(0)  is  still  estimated  with  the  full  degrees  of  freedom,  but  R^t) 
for  t  0  is  estimated  with  fewer  degrees  of  freedom.  For  t  near  L2,  the  loss  in 
degrees  of  freedom  in  estimation  of  R,(t)  is  significant,  and  the  variance  increases. 
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Figure  16.  Normalized  Quality-Ratio  for  Overlapped  Rectangular 
Temporal  Weighting  and  Lag  Reshaping  to  Cl 
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q  =  .5  50  Percent  Overlap 

For  P  =  l,  there  is  no  overlap  and  conditions  are  identical  to  those  described 
above.  For  P  =  2,  the  sudden  increase  in  variance  can  be  explained  as  follows:  from 
(16),  the  first-stage  estimate  is 


A 

Rj(t) 


-?[/ 


dt  y0(t)  Xo(t  " 


yjCt) 


In  particular, 


^(0) 


-[/< 


|yQ(t) 


Ixjtt) 


(145) 


where  the  overall  weighting  O(t)  of  x2(t)  is  depicted  in  figure  17.  As  shown,  the 
overall  weighting  is  very  uneven,  causing  loss  in  stability.  As  P  increases  above  2, 
the  uneveness  of  the  overall  weighting  (for  q  =  .5)  occurs  only  towards  the  edges  of 
the  (0,  T)  interval,  yielding  a  decrease  of  variance,  since  more  data  points  tend  to  get 
the  same  overall  weighting,  insofar  as  their  effect  upon  the  estimation  of  R,(t)  is 
concerned.  However,  at  the  same  time,  the  effect  of  fewer  degrees  of  freedom  in 
estimation  of  R,(r),  for  r  values  near  L2,  becomes  more  pronounced  as  P  increases 
and  L,  decreases;  this  is  true  even  for  the  50  percent  overlap  case  being  considered 
here.  Eventually,  this  effect  dominates,  and  the  variance  increases  with  P. 
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q  =  .75  25  Percent  Overlap 

For  any  value  of  P  >  1,  the  overall  weighting  O(t)  is  very  jagged  (as  above),  and 
the  jaggedness  does  not  decrease  or  concentrate  near  the  edges  as  P  increases.  This  is 
true  for  any  overlap  greater  than  0  and  less  than  50  percent. 

In  summary,  for  rectangular  temporal  weighting,  the  smallest  values  for  the 
normalized  quality  ratio  are  realized  by  choosing  q  =  1,  no  overlap,  and  making  L| 
several  times  larger  than  L2.  This  conclusion  about  the  ratio  L^/L,  is  consistent  with 
those  reached  earlier. 

Hanning  Temporal  Weighting 


The  temporal  weighting  and  associated  correlation  for  this  case  are  given  by 
(116),  (C-9),  and  (C-10): 


Evaluation  of  third-order  correlation  $3(t,  p)  in  (108)  is  rather  tedious;  the  end 
result  is  given  in  (F-l)  and  (F-2).  The  procedure  and  program  for  the  evaluation  of 
the  normalized  quality  ratio  is  given  in  appendix  F. 


The  normalized  quality  ratio  for  Hanning  lag  weighting  and  BeT  =  100  is  plotted 
in  figure  18  for  several  values  of  the  shift  fraction  q.  When  q=  1,  no  overlap,  the 
small  values  of  the  Hanning  temporal  weighting  at  its  edges  cannot  be  compensated 
for,  by  any  choice  of  L2/L,,  and  the  variance  remains  at  approximately  twice  the 
ideal  value.  For  50  percent  overlap  of  the  Hanning  temporal  weighting,  q  =  1/2,  the 
situation  is  markedly  improved,  there  being  a  value,  L2/L,  =  .4,  at  which  the  excess 
variance  is  only  8  percent  above  ideal;  this  is  reminiscent  of  the  variance  ratio  for 
the  case  of  no  quadratic  smoothing  in  reference  9,  tables  5-7.  When  q  is  decreased  to 
3/8  or  1/4  (62.5  and  75  percent  overlap,  respectively),  virtually  the  idea!  variance 
reduction  can  be  achieved  by  choosing  L2  =  L,. 

In  figure  19,  the  shift  fraction  q  is  kept  at  3/8,  while  BeT  is  taken  at  both  100  and 
32.  The  smaller  value  of  BeT  leads  to  a  slightly  larger  loss  in  performance  because  of 
more  significant  edge  effects.  However,  even  so,  the  normalized  quality  ratio  does 
reach  a  very  desirable  level  only  4  percent  above  ideal  when  L2  ~  2L,. 
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Figure  19.  Normalized  Quality-Ratio  for  Hanning  Temporal  Weighting 
and  Hanning  Lag  Weighting;  q  =  3/8 
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Summary 

The  possibilities  and  performance  of  a  generalized  spectral  analysis  technique 
employing  temporal  and  lag  weighting  have  been  investigated  in  terms  of  the  mean 
and  variance  of  the  spectral  estimate.  The  only  assumption  required  about  the 
process  under  analysis,  in  so  far  as  the  mean  is  concerned,  is  that  it  be  second-order 
stationary  over  the  observation  interval.  We  then  were  able  to  extract  a  simple 
expression  for  the  effective  window  involving  the  temporal  and  lag  windows. 

The  possibility  of  doing  lag  reshaping  to  ach^ve  desirable  effective  windows  was 
considered  in  detail  and  found  to  be  reasonable  for  a  wide  variety  of  windows  with 
good  side  lobe  behavior  and  decay  rates.  In  particular,  if  rectangular  temporal 
weighting  is  employed,  its  inherent  poor  side  lobe  structure  can  be  corrected  via 
proper  lag  weighting,  in  so  far  as  the  effective  window  is  concerned.  Strictly 
speaking,  the  closest  side  lobe  cannot  be  eliminated;  however,  all  the  other  side 
lobes  can  be  suppressed. 

The  effect  of  temporal  and  lag  weighting  on  the  variance  of  the  spectral  estimate 
was  evaluated  and  compared  with  the  ideal  value  for  large  BeT.  For  rectangular 
temporal  weighting,  it  was  found  that  small  values  of  L2/Lj  and  no  overlap  led  to 
values  of  the  normalized  quality  ratio  virtually  equal  to  the  best 
attainable  by  any  spectral  analysis  technique.  The  comparison  is  made  under  the 
constraint  that  the  effective  frequency  resolution  Be  is  maintained  the  same  for  all 
techniques  under  consideration.  On  the  other  hand,  if  Hanning  temporal  weighting 
is  employed,  overlapping  must  be  used  for  maximum  variance  reduction  and  the 
length  ratio  L2/L,  ought  to  be  of  the  order  of  unity. 

Since  Hanning  temporal  weighting  requires  multiplication  of  each  and  every  data 
segment  (P  pieces)  and  significant  overlap  (~  50  percent),  whereas  rectangular 
temporal  weighting  requires  no  multiplication  and  no  overlap,  the  latter  approach  is 
a  strong  candidate  for  spectral  analysis,  particularly  since  excellent  effective 
windows  (low  side  lobes  and  rapid  decay)  and  virtually  ideal  variance  reduction  can 
be  achieved  by  proper  lag  weighting  and  choice  of  ratio  L2/Lr  Investigation  of 
other  cases  than  those  evaluated  here  can  be  achieved  by  appropriate  modification 
of  the  program  in  appendix  F.  A  major  analytical  task  will  be  the  evaluation  of  the 
third-order  correlation  (108),  if  temporal  weighting  Wj(t)  is  taken  other  than  rec¬ 
tangular  or  Hanning. 
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Appendix  A 

Comparison  of  Two  Bandwidth  Measures 


The  effective  bandwidth  of  narrowband  window  W  (f)  was  defined  in  (1)  as 


B  = 

e  ( 


df  WQ(f) 


[ df  W2(f) 

J  o 


The  half-power  bandwidth,  Bh,  is  defined  as  the  solution  of 


(fo  -  l  Bh)  1 


W 


where  it  is  assumed  that  window  Wo(0  is  real,  even  about  f0,  and  peaked  at  f0.  We 
let  We(f)  =  W0(f  +  f0);  thus  We(f)  is  even  about  f  =  0. 

The  inverse  Fourier  transform  of  lowpass  window  We(f)  is  called  the  weighting 


w  (t)  =  f df  exp(i2Trfx)  W  (f) 

c  J  0 


We  consider  here  the  class  of  weightings  given  by 


w  (t)  =  cos 

e  k>0  K 


^irkt/L^  for  |t|  < 


and  zero  otherwise,  where  {ak}  are  real  and  non-negative.  This  class  includes 
rectangular,  Hanning,  Hamming,  Blackman,  and  the  optimal  windows  of  Nuttall, 
reference  19.  The  Fourier  transform  of  (A-4)  yields  lowpass  window 

2Le  ^  (-1)K  ak 

W  (f)=  —  2L  f  sin(2rrL  f)  £  - 7  - -  -•  ; 

e  "  e  e  k>0  (2L  f)2  -  k2 

e 


Mt)  ■ 


2aQ  for  n  =  0 


a  i  |  for  n  ^  0 


A  table  of  bandwidths  Be  and  Bh  and  their  ratio  is  given  below  for  the  window  in 
(A-5).  Although  these  bandwidths  vary  significantly  for  the  different  weightings, 
their  ratio  is  much  more  stable.  In  fact,  for  the  last  four  weightings  listed,  the  ratio 
is  constant  within  ±  1  percent.  The  weightings  listed  under  C5,  C3,  Cl  are  those 
given  in  reference  19,  figures  10, 11, 12;  the  notation  means 


A-l 
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C5:  continuous  fifth  derivative  of  weighting 
C3:  continuous  third  derivative  of  weighting 
Cl:  continuous  first  derivative  of  weighting 


Table  A-l.  Bandwidths  for  Various  Weightings 


Weighting 

BeLe 

Rectangular 

0.5000 

Hanning 

1.3333 

Hamming 

1.2614 

Blackman 

1.6415 

C5 

2.2165 

C3 

2.0478 

Cl 

1.9544 

BhLe  Be/Bh 


0.6034 

0.8287 

1.0000 

1.3333 

0.9109 

1.3848 

1.1494 

1.4281 

1.5371 

1.4420 

1.4139 

1.4483 

1.3444 

1.4537 
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Appendix  B 

Some  Lag  Weighting  and  Lag  Windowing  Considerations 
For  Discrete-Time  Processing 

It  is  convenient  here  to  define  an  equispaced  unit-impulse  train  by  the  notation 


6  (b)  =  y  6  (b  -  na)  , 
av  n  (B-l) 

where  the  summation  on  n  extends  over  ±  00 . 

For  discrete  time  sampling  at  spacing  A,,  it  has  been  observed  in  (30)  that  6,(f) 
has  period  1/A,  in  f.  Therefore  lag  window  W2(f)  could  be  confined  to  |f|  <  (2A,)-' 
with  no  loss  in  generality,  in  so  far  as  its  effects  on  62(f)  by  means  of  (18)  are 
concerned.  In  fact,  for  a  general  lag  window  W2(f)  specified  arbitrarily,  the 
equivalent  band  limited  lag  window  is 


W2(f)  =  rect(Atf) 


w2(f) 


9  6  t  (£) 


where  we  will  utilize  definitions 


rect (x) 


1  for  j  x |  <  1/2 
0  otherwise 


sine (x) 


sin(irx) 

irx 


(B-2) 


(B-3) 


and  where  ©■  denotes  convolution.  That  is,  W2(f)  is  aliased  into  the  band 
|  f |  <  (2A,)-1,  and  only  this  band-limited  portion  is  retained  for  $ 2(f ). 

A  way  to  demonstrate  this  mathematically  is  to  note  that  the  only  values  of  lag 
weighting  w2(t)  that  can  affect  &2(t)  are  the  samples 


w2(qAt)  =  Jdf  exp(i2nfqAt)  W2(f) 


(B-4) 


The  band-limited  lag  weighting  function  that  passes  through  all  these  specified 
values,  for  all  q,  is 


w2(t) 


Ew2(q  V 

q 


‘(B-5) 


with  corresponding  Fourier  transform 
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W2(f)  =  A  rect(Atf)  £ 

q 


w2(qAt)  exp(-i2irfqAt) 


(B-6) 


=  A  rect(Atf) 


^2  exp(- i2irfqAt)  J  du  exp(i2TruqAt)  W2(u) 
q 


j  du  W2 (u)  At  £  exp^-i2Tr(f  -  u)qAt 


j  du  W2 (u)  6  j  (f  -  u) 


=  rect(Atf) 


Relation  (B-4)  indicates  how  an  arbitrarily  specified  W2(f)  fixes  the  lag  weights  at 
the  sample  points.  For  the  reverse  problem,  where  sampled  lag  weights  {w2(qA;)} 
are  specified  for  all  q,  relation  (B-6)  gives  an  equivalent  lag  window,  in  particular 
the  band-limited  spectral  window,  which  results  in  the  same  estimates  R,(t)  and 
62(f).  Notice  that  W2(f)  is  not  uniquely  specified  by  samples  {w2(qA,)};  however, 
the  band-limited  W2(f),  which  realizes  weights  |w2(qA,)}  for  all  q,  is  unique  and  is 
given  by  (B-6). 


w2(f) 


9 


6  j  (f) 


(B-7) 


=  rect(Atf) 


=  rectfA^f) 


As  a  special  case  of  the  above,  consider  discrete  frequency  smoothing  with 
frequency  spacing  Af  =  (MA,)-1,  where  MA,  is  of  the  order  of  2L , ;  i.e.,  from  (38), 


W2^  =  MA,  ^ 


2n 


(B-8) 


where  we  set  W2n  =  0  for  |n|  >  M/2  without  loss  of  generality,  in  accordance  with 
the  observation  above  (B-2).  Then  lag  weights  (B-4),  given  now  by 


W2(<*V  =  MK7 


n  I  <M/2 


^2n 


(B-9) 


will  equivalently  accomplish  the  same  purpose.  This  last  relation  can  be 
accomplished  by  an  M-point  FFT,  where  W2  ±M/2  receive  the  same  complex 
exponential  weighting  in  (B-9). 

It  should  be  noted  that  the  discrete  function  w2(qAt)  in  (B-9)  has  period  M  in  q; 
this  means  that  w2(qA,)  in  (B-9)  will  increase  in  magnitude  for  M/2  <  q  <  M.  If  ft  , 
is  nonzero  for  |q|  >  M/2,  this  lag  weighting  may  cause  a  problem.  One  guaranteed 
way  to  avoid  the  problem  is  to  choose  M/2  larger  than  the  nonzero  extent,  N,,  of 
ft,q.  Physically,  this  means  that  the  frequency  spacing  Af  =  (MA,)-1,  used  in 
frequency  smoothing  (B-8),  must  be  small  enough  so  as  not  to  miss  any  information 
in  6,(0-  Coarse  frequency  spacing  gives  spurious  results  for  62(f).  (It  will  also  yield 
poor  effective  windows.) 


B-2 
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Since  from  (24),  k  ,q  =  0  for  q  ^  N,,  where  N,  is  the  number  of  time  samples  per 
segment,  only  a  finite  number  of  the  general  weights  (w2(qA,)}  in  (B-4)  affect  R2q. 
Thus  in  example  (B-8),  although  (B-9)  has  period  M  in  q,  only  the  values  for 
|q|  <  N,  are  relevant  to  the  effect  on  &2q;  more  generally,  the  values  yielded  by  (B-4) 
for  a  general  W2(f)  are  relevant  only  for  |q|  <  N,,  and  only  these  need  be  evaluated 
and  retained  if  we  choose  to  process  via  the  lag  domain. 

Now  let  us  consider  the  reverse  problem,  where  lag  weights  {w2(qA,)}  are 
specified  for  all  q,  and  we  wish  to  determine  some  allowable  lag  windows  W,(f)  that 
will  realize  the  same  estimates  ft2(r)  and  <J2(f),  but  which  take  advantage  of  the  fact 
that  only  w2(qA,)  for  |q|  <  N,  must  be  realized.  One  obvious  candidate  is  the  band- 
limited  lag  weighting  version 

w  (t)  =  £  w  (qA  )  sincnr  '  °l)  • 

|q|<N1  \At  / 

notice  the  limitation  on  q  employed.  The  corresponding  lag  window  is 

W2(f)  =  A  rect(A  f)  5Z  w2(qAt) 

Iq|<Nj 

A  second  candidate  is 

S2(t)  =  w2(t)  ®  5^  (t) 

provided  that  M  >  2N,-1;  this  provision  guarantees  non-overlap  of  the 
displacements  of  \V2(t).  Then 

MS'  '  £  W2n  S(f  -  MS')  (B.I3) 

z  MA  n  X  7 

where 

W2n  =  mZT  ^2  (m ?-)“£  rect(£),  £,  W2(qAt3  exp(-i2wnq/M) 

t  \  t/  ' |q| <Nj 

(B- 14) 

from  (B-ll).  Notice  that  (B-13)  has  the  form  of  discrete  frequency  smoothing  in 
(B-8);  (B-14)  gives  the  area  of  each  impulse  needed  in  (B-13).  Also  notice  from 
(B-14)  that  all  these  areas  are  zero  for  |n|  >  M/2;  thus  we  have  a  finite  sequence  to 
apply  in  the  frequency  domain,  which  is  equivalent  to  a  specified  finite  set  of  lag 
weights. 

Equations  (B-14)  and  (B-9)  are  complementary  to  each  other.  In  fact,  we  can 
derive  (B-9)  from  (B-14)  as  follows:  from  (B-14), 

W-  =  ^  rectfrr)  y*  w  exp(-i2nnq/M)  , 

2n  M  {Ml  Lf  q  ^ 


(B-10) 

(B-ll) 

(B-12) 
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where 


w  = 

q 


Then 


|n|<M/2 


where 


•  CM) 


(  w2(qAt) 

for  [ q |  <  Nj 

(  ° 

otherwise 

-  £ 

exp(i27rnp/M)  * 

|n  1  <M/2 

*  £  "q 

1  o 

m  £  rect 

M  1  n  |  <M/  2 

T  q 

*  £  w 

rW  -  w  ®  i 

t-j  q 

q 

p-q  p 

for  p 

=  0,  +  M,  _+  2M, 

0  otherwise 

Now  if  M  >  2N.-1  (as  assumed  above  (B-13)),  then 


£  W»  exp(i2irnp/M)  =  w?(qA  )  for  |q|  <  N. 
|n|<M/2  zn  1  z  1 


This  is  (B-9)  for  |q|  <  N ,,  which  is  the  only  range  that  affects  6 


2q‘ 


B-4 


(B-16) 

exp(-i2imq/M) 

-  q)/M) 

(B-17) 

(B-18) 

(B- 19) 
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Appendix  C 

Correlation  $,(t)  of  a  General  Class  of 
Temporal  Weightings 

The  class  of  temporal  weightings  of  interest  here  is  given  by  a  sum  of  complex 
exponentials: 

w.(t)  =  53  ch,  exp(i2Trkt/L1)  for  |tj  <  L./2  , 

k  r-1) 

and  zero  otherwise.  We  assume  that  the  coefficient  sequence  has  conjugate  sym¬ 
metry 

“  k  =  ak  ^or  a11  k  * 


then  w.(t)  is  real,  and  it  follows  that  the  (aperiodic)  correlation 


>2  Ct)  =  J dt  Wj(t)  Wj(t  -  t) 


is  also  real,  in  addition  to  being  even  about  t  =  0. 
Substitution  of  (C-l)  in  (C-3)  yields 


L2/2 

^(t)  =  £  «*  J  exp^i2irkt/Lj  -  i2irm(t  -  t)/LjJ 


t-Lj/2 


for  0  <_  t  <_  Lj  .  (C-4) 


This  can  be  evaluated  and  then  extended  to  r  <  0  by  the  use  of  the  even  character  of 
4|(t);  there  follows 

^1  (T)  =  *  M)  |uq  +  2  53  l<*k|2  cos(2Ttki/L1)j 


2L1  v  (-l)m-k 


-  ~r  £ 


^ak  a*^  sin(2Ttm|x|/L1) 

for  | t |  <  L1  ,  (C-5) 


and  zero  otherwise.  This  is  the  general  result  for  the  correlation  of  weighting  (C-l ). 
We  now  specialize  (C-5)  to  the  case  of  real  symmetric  coefficients  in  (C-l): 

°k  real  ’  a-k  =  °k  * 


C-l 
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For  all  coefficients  zero  except  for  cr0,  we  have 

Wj(t)  =  aQ  for  | t |  <  Lj/2  , 


(C- 7) 


«  (Lj  *  |t|)  Oq  for  |t|  <  Lj  . 

For  all  coefficients  zero  except  for  ag,  at,  we  have 

Wjft)  =  «q  +  2a1  cos(27rt/Lj)  for  |t|  <  Lj  /2  , 

=  (Lx  -  |  t  | )  [4  +  2a^  cos^inr/l^)} 

L1 

+  —  a1(2a0  "  aj)  sin(2n |t |/Lj)  for  | t 1  <  Lj 


(C-8) 


(C-9) 


(C-10) 


For  all  coefficients  zero  except  for  ao,  a,,  o2,  we  have 

Wj(t)  =  aQ  +  2a L  cos(2irt/L *  2a 2  C0s(4»t/Lj)  for  Jt  J  <  Lj/2  , 

(C-ll) 


q  +  2a^  cos(2ht/Lj)  +  2a^  cos(4irT/Lj)J 

q  -  3at  -  4a2)  sin(2it |t \/L^) 

16aj  +  3a2)  sin(4ir  1 1  l/L^J 

for  | t |  <  Lj  •  (C-12) 


For  all  coefficients  zero  except  for  a0,  o,,  a2,  ay  we  have 

Wj(t)  =  +  20^  cos^nt/Lj)  +  2a2  cos(4iit/Lj)  ♦  2a3  cos(6nt/Lj) 

for  |t|  <  Lj/2  ,  (C- 13) 


C-2 


^(t)  =  (Lx  -  I  T  I  ) 

+  2a^  cos (2-itt/Lj)  +  2a^  cos(4itt/Lj)  +  2a^  cos (6ttt/L1)J 

l*  r 

+  ISajCnap  -  6otj  -  8a2  +  3a3)  sin(2ir  |  t  | /l^) 

-  a2(30aQ  -  80a^  +  15a2  +  48a3)  sin(4ir | t |/Lj) 

+  a3(20a0  -  450^  +  72a2  -  10a3)  sin (6ir | x  | /Lj)J 


for  lx  I  <  L 


1 


(C-14) 


This  last  case  includes  all  the  weightings  considered  in  reference  19,  with  the 
identification  of  coefficients  as 


a0  a0 


ak  =  yak  for  k  >  0  . 


(C-15) 


Then  we  can  express  the  temporal  weighting  as 

w- (t)  =  £  av  cos(2irkt/L1)  for  |t|  <  L 1  / 2 
k>.0  K 


(C-16) 
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C-3/C-4 
Reverse  Blank 
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Appendix  D 

Derivation  of  Variance  of  Spectral  Estimate 

Our  starting  point  is  (14).  The  integral  on  t  is  a  Fourier  transform  of  the  product 
in  (13),  and  can  therefore  be  expressed  as  a  convolution: 

Yp(f)  *  X(f)  ®  j^W^f)  exp^-i2irf^-y-  +  psjj 

*  J  du  X(f  -  u)  Wj(u)  exp +  psj^  , 


(D-l) 


where  we  used  (46)  and  defined 

X(f)  =  J dt  x(t)  exp(-i2irft) . 


(D-2) 


Although  the  relations  to  follow  could  be  derived  in  the  time  domain,  it  is  more 
convenient  to  develop  them  in  the  frequency  domain  because  of  the  frequent  and 
useful  occurrence  of  delta  functions. 

X(f)  is  complex  Gaussian  for  all  f,  since  (D-2)  is  a  linear  transformation  and  we 
have  assumed  x(t)  to  be  a  complex  Gaussian  process,  for  the  variance  calculation  to 
follow.  Furthermore  covariance 

Av{X(fj)  X*(f2)}  =  JJ dtj  dt2  x(tj)  x*(t2)  exp(-i2irf jtj  +  i2irf 2t^) 
=  ffdt j  dt2  R(tj  -  t2)  exp(-i2nf1tj  +  i2irf2t2) 

=  JJ  du  dt2  R(u)  exp(-i2irf  1  (u  ♦  t2)  +  i2irf2t2) 


=  Gffj)  d(fj  -  f2) 


CD-3) 


upon  use  of  (43)  and  (45).  When  x(t)  is  a  single-sided  (analytic)  complex  process, 
there  then  follows  for  the  fourth-order  average,  which  will  be  needed  later 
(reference  23), 


Av { X ( f : )  X*(f2)  X(f3)  X*(f4)> 

=  G(f j)  G(£3)[fi(f1  -  f2)  «(f3  -  f4)  4  «(fj  -  f4)  <S(f 2  -  f3)] 

(D-4) 


L 
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(When  x(t)  is  a  real  process,  there  would  be  a  third  term  in  addition.  However,  even 
then,  this  term  contributes  only  near  f  =  0;  see  reference  9,  equations  (A-4)  -  (A-8).) 

Upon  use  of  (D-l),  (14)  becomes* 

A  i  P”1  rr 

Gt(f)  =  p  E  II  du  dv  x(f  -  u>  x*(f  -  v) 

p=0  JJ 


x(u)  W* (v)  expl -i27i(u  -  v)l -j-  +  pS 


Then  the  average  of  the  product  of  the  first-stage  spectral  estimates  is 

P-1 

AvfGjffj)  G^  (f  2  D }  =  ^2  £  (  ,  ////  du  dv  du  dv 

xlT;  -  u)  X*(ft  -  v)  X(f2  -  u)  X*(f2  -  Z)  W1( U)  Wx*(v)  Wx(u)  W*(v) 


•  expl  -i27r(u  -  v) 


&  *  m)' 


expl -i2ir(y  -  v)\—  +  qS 


Reference  to  (D-4)  enables  us  to  express  the  fourth-order  average  as 

G (f  -  u)  G(£2  -  y) 

(u  -  v)  <5(y  -  v)  +  Sffj  -  f2  +  v  -  q)  6(fj  -  f2  +  y  -  v)]  .  (D-7) 

Use  of  the  first  term  of  (D-7)  in  (D-6)  yields 

i  rr  2  2 

X  £  JJdu  dy  G(f T  -  u)  G(f  j  -  V)  IWjfu)!  JWjCw)  | 

P  p,q=0  JJ 

=  J  du  GCf^  -  u)  |WX  (u)  l2  •  J  dy  G(f  2  -  u)  |  W^,  ( yiD  | 2 

=  Av{5j  (f  j)  }  •  Av{^1  (£2)  >  ,  (D  8) 

where  we  employed  (45)  in  the  last  line.  Moving  this  term  to  the  left  side  of  (D-6), 
and  using  the  second  term  of  (D-7),  we  obtain,  for  the  covariance  of  the  first-stage 
spectral  estimates. 


•For  more  generality,  we  allow  temporal  weight  w(  and  window  W(  to  be  complex  in  this  appendix. 


D-2 
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c°v(G  (f  ),  GjCfJi  =  X  2:  ffdudv  G(fj  -u)  G(fj  -v) 

11  1  4  P p, q=0  1  J 

Wj (u)  W* (v)  Wj(v  +  f2  -  fj)  W|(u  +  f2  -  fj)  exp(-i2Tr(u  -  v)  (p  -  q)S) 

1  P -1  f  f 

=  -7  £  J  J  da  d6  G M  G(6)  W1(f1  -  a)  W*(f1  -  6) 

P  p,q=0 

*  W1(f2  -  0)  W*(f2  -  a)  exp(-i2jr(0  -  a)  (p  -  q)S) 

,  P-1  1  ,  2 

=  -j  E  J  da  G(a)  WjCfj  -  a)  W*  (f 2  -  a)  exp(i2™(p  -  q)S) 

P  p  ,q=0  I 

J  da  G(a)  Wj(f^  -  a)  W*(£^  -  a)  exp(i27iapS) 


p=l-P 


=// da  d0  G(a)G(6)  WjCfj  -  a)  W*^  -  0)  Wj(f2  -  6)  W*(f2  -  a) 


[sin(TrPS(a  -  0) ) 

P  sin(7rS(a  -  0)}J 


Here  we  have  used  the  identities 


1 

“2 


P-1 

E 

p,q=0 


exp(i27r(p 


1  P'1 

q)“)  =  f  E 

P=l-P 


(D-9) 


exp(i27Tpu) 


sin  (ttPu) 
P  sin(7TU 


;  Qp(u) 


(D-10) 


For  f,  =  f2,  (D-9)  checks  with  reference  9,  equation  (A-9);  more  generally,  it  is 
equation  D-2  of  reference  9.  We  observe  that  if  |f2  -  f,|  is  greater  than  the  effective 
bandwidth  of  temporal  window  W,,  (D-9)  will  be  small  since  W^f,  -  0)  and 
w;(f2  -  0 )  will  then  not  overlap  significantly  on  the  0-scale.  Also  notice  that 
spectrum  G  is  still  left  under  the  integral  sign;  i.e.,  there  are  no  assumptions  yet  on 
the  character  of  the  spectrum. 

We  are  now  prepared  to  consider  the  second-stage  spectral  estimate  62(f)  as  given 
by  (18): 

G2(f)  =  fdf1  GjCfj)  W2(f  -  fj)  .  (D-Il) 


Then 


D-3 
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S2(f)  -  Av{G2(f)}  =  /df^Vfj)  -  AvjGjCfpjJ  W2(f  -  fj)  , 


(D-12) 


and  therefore  (recall  that  W2  is  real)  the  variance  of  the  second-stage  spectral 
estimate  is 

Varjc2(f)j  =  ffd£1  df2  CovJCjCfj).  W2(f  -  «2(f  -  ff2) 

da  d 6  G(a)  G(6)  exp(i2?r(a  -  B)pS)|y(f  -  a,  f  -  8) 

=  f  j da  dg  G(a)  G(6)'  |y(f  -  a,  f  -  B)  |2  Qp(S(a  -  B))  ,  (D.,3) 

where  we  used  (D-9)  and  (D-10),  interchanged  integrals,  and  defined  window 
convolution  function 


y(x,y)  =  J du  W2(u)  Wjfx  -  u)  W*(y  -  u) 


(D-14) 


Relation  (D-13)  is  exact;  it  makes  no  presumption  about  the  relative  widths  of  the 
spectrum  G  and  the  windows  WJ#  W2.  The  compact  expression  (D-13)  involves  the 
windows  Wp  W2  through  the  convolution  function  y,  and  involves  the  shift  S  and 
number  of  pieces  P  through  the  periodic  function  QP  defined  in  (D-10). 

The  window  convolution  function  y  in  (D-14)  realizes  its  peak  value  at  x  =  0,  y  =  0, 
and  is  rather  small  everywhere  else,  since  the  windows  are  virtually  unimodal  and 
rather  narrow.  In  fact,  a  special  case  is  the  diagonal  slice 

y(x,x)  =  f du  W2(u)|W1(x  -  u)|2  =  |wi(x)|2  9  W2(x)  =  Wg(x)  , 

(D-15) 

by  reference  to  (51).  Generally,  y(x,  y)  is  substantially  nonzero  only  in  the  region 
Be,  Be  at  the  origin  of  x,  y  space. 

We  now  employ  the  assumption  discussed  in  connection  with  approximation 
(103),  namely,  that  spectrum  G  is  relatively  constant  in  the  band  of  width  Be  about 
the  frequency  f  of  interest.  Then  the  major  contribution  to  the  variance,  (D-13), 
comes  from  the  region  near  a  =  p  =  f  in  a,  p  space.  There  follows  the  approximation 
for  the  variance  of  the  second-stage  spectral  estimate, 


D-4 
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Var{GA2(f))  .  G2(f)  p  (l  -  M) 

•  J J  da  d6  exp(i2u(ai  -  8)pS)|y(f  -  a,  f  -  B)  | 2 
=  c2(f)  i  £  p(-¥)//  dx  dy  exp(i27t(y  -  x)pS)  | y (x ,  y) 


=  G2(f)  j  j  dx  dy  |y(x,  y) |2  Q  CS (y  -  x))  , 


(D-16) 


where  we  used  (D-10)  again. 


We  now  simplify  the  double  integral  in  (D-16);  from  the  second  line  of  (D-16)  and 
from  (D-14), 

D  =  J  J  dx  dy  exp(i2ir(y  -  x)pS)  J  J  du  W^fu)  W^(x  -  u)  W*(y  -  u)  J2 

=  J  J  du  dv  W2(u)  W2(v)Jjdx  exp(-i2irXpS)  Wj(x  -  u)  W*(x  -  v)|2 

=  //  du  dv  W2(u)  W2(v)  IxjCpS,  u  -  v)  |2  >  (D-17) 

where  the  complex  ambiguity  function  of  the  temporal  weight  and  window  is 
defined  as 

xi(t,  v)  =  J df  exp(i2TrfT)  W^f  +  ^  W*^f  - 

=  / dt  exp(-i2™t)  Wl(t  +  J)  w*(t  -  J-)  .  (D  lg) 

Now  let  v  =  u  -  v  in  (D-17),  and  obtain  a  single  integral  for  D: 

D  =  J  dv | xx (pS,  v)|2  J  du  W2(u)  W2(u  -  v) 

=  /  dvIxjfpS,  v)|2  x2(0,  v)  ,  <D-19) 

in  terms  of  the  ambiguity  function  y2  of  the  lag  window  and  weight.  Substitution  of 
(D-19)  in  the  second  line  of  (D-16)  yields  for  the  variance, 

Var{G2(f)}  =  G2(f)  ~  E  (l  -  fdv|Xl(pS,  v)|2  x2(0,  v)  . 

P=1-pV  (D-20) 

An  alternative,  and  perhaps  more  useful  form,  to  (D-20)  is  attained  as  follows; 
the  integral  in  (D-20)  is  expressible  (by  use  of  the  definition  (D-18)  of  the  ambiguity 
function)  as 


D-5 
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Jdvfx^pS,  v)T  /  dt  exp(-i2irvt)  w^(t) 

=  Jdt  w22(t)  J  dv  exp(-i2irvt)  |  Xj  (pS,  v) 


(D-21) 


Now  we  have  the  general  result  that 

fdv  exp(-i27Tvt)  | x1  Cx ,  v)  |  2 
=  j  dv  exp(-i27rvt)  d*2 


exp(-i2TTv  (t^ 


t2)) 


•  fdti  ui('i  *  7)  ";(ii  -  f) »;('  **1*7)  "if  *  'i 
/ du  “if  *  i~ts)  “if  4  tLrJ~)  “if  *  LiJ-)  “if  * 


-T  +  t 


=  4>3Ct.  t)  .  (D-22) 

This  is  a  third-order  correlation  of  temporal  weighting  w,;  see  the  footnote  to  (9). 
Thus  (D-2 1 )  becomes 


Jdt  w^Ct)  <f>3(t,  pS) 


(D-23) 


and  the  variance  in  (D-20)  becomes 

Var(G*2(f)J  -  C2(f)  i  f  -  M)/dt  w2(t)  *,(t,  pS) 

This  “weighting-domain”  version  of  the  variance  is  very  useful  if  third-order 
correlation  can  be  evaluated  in  closed  form. 
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Appendix  E 


Computational  Considerations  for  Non-Overlapping  Segments 


The  curves  in  figure  15  for  non-overlapping  time  segments  are  drawn  over  a 
continuum  of  values  of  L2/Lt.  However,  if  we  were  given  a  value  of  TBC,  all  values 
of  L2/L,  may  not  be  allowed.  To  develop  this  point,  suppose  that  we  pick  an  integer 
value  for  the  number  of  pieces,  P,  and  solve  for  L,Be  =  TBe/P  according  to  (129). 
From  figures  5  and  6,  this  dictates  the  value  for  L2Be,  and  hence  a  discrete  value  for 
L2/L,  is  specified  for  each  value  of  P.  The  number  of  pieces,  P,  can  range  from  1  to 

Pmax.  where 

B  T 

B  L  (min)  =  =  2c{V"}  ’  Pmax  =  2c{4>1}  Be  T  *  (E-l) 

max  1 


from  figure  5.  For  BeT  »  1,  we  have,  for 


1>  BeLl  =  Be  T’  BeL2  *  2c  {w“}  *  *  2^Tb^T 


<<  1 


(E-2) 


and  for 


P 


=  P 


max 


BL  = 
e  1 


2c{<(,1\  * 


B  L,  >  B  T, 
e  2  —  e 


T  >>  1  . 

(E-3) 


Thus  a  very  wide  range  of  discrete  values  of  L2/L,  is  allowed  when  BeT  »  1 . 

The  problem  with  this  approach  is  that  when  BeL,  is  calculated,  BeL2  must  be 
solved  for  from  the  integral  relation  (52)  (or  approximately  from  figure  6).  This 
tedious  procedure  can  be  circumvented  by  specifying  L2/Lj  instead;  if  desired,  we 
could  then  use  (52,  to  determine  BeL)(  and  solve  for  P  =  BeT/(BeL,).  However,  P 
will  not  necessarily  turn  out  to  be  an  integer  for  a  given  fixed  T;  thus  only  a  set  of 
discrete  values  of  L2/L,  are  strictly  legal.  But  if  NQRPin  (130)  does  not  vary 
radically  with  L2/L,,  this  is  not  a  significant  limitation.  And  since  it  is  simpler,  we 
adopt  it.  In  Appendix  F,  we  cannot  avoid  the  calculation  of  BeL2  from  a  given  BeL , . 


E-l /E-2 
Reverse  Blank 


Appendix  F 


Computational  Considerations  for  Overlapping  With 
Hanning  Temporal  Weighting 

The  temporal  weighting  is  given  by  (146).  Substitution  in  (108)  and  evaluation  of 
the  integral  yields 

4  k-1 

$3(t,  u)  =  g^r—  Vq(w  -  a  -  B)  +  £  (_1)  vk  sin(k(a  +  B)) 

L1  k=l 

for  a  +  B  <.  it  ,  (F-l) 


where 


a  =  f|x|  ,  6=  f|u! 

L1  L1 


=  cos(a)  =  coshf-j  »  =  c°s(S)  =  cosf~ 


*22  2  2  4  4 

V0  =  8  "  Ca  -  CB  +  4Ca  C8  +  C«  +  C8 


V,  =  C  C.(4C2  +  4C2  -  1 
1  a  3  \  a  3 


v?  =  -  j  +  yfc2  +  c«)  +  c2  cft  , 

2  4  2  \  a  B/  a  6 


V2  "  3Ca  CB'  V4  =  32  ‘  (F-2) 

The  procedure  for  the  evaluation  of  the  normalized  quality  ratio  follows.  We 
specify  a  value  for  BeT  and  select  a  temporal  weighting  w,(t)  and  a  lag  weighting 
w2(t).  We  then  evaluate  shape  factor  c{+, }  from  (64)  or  table  2,  and  select  a  shift 
fraction  q  -  q  {w,}  according  to  (139).  We  then  solve  (141)  for  Pmax,  and  allow  P  to 
take  integer  values  in  the  range  1  <  P  <  Pmax.  B^  can  then  be  evaluated  from 
(140)  as 

B  T 

BeLl  =  1  ♦  (P  -  l)q 


4  *  f ufif  >jlur.»  P  U  c  fhcn  S4?|\c  (?2)  in  the  form 


[•t  w'n.jO  fi  x)  , 


?  !  ♦  <1  xmn  independent  ot 
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(F-5) 


where  we  let  t  =  L,x  and  used  (139).  The  quantity  w2(L,x)  is  a  function  of  only 
L-./L,,  while  the  remaining  quantity  in  the  integrand  of  (F-5)  is  independent  of  L, 
and  L2.  Finally,  we  multiply  (F-5)  by  BeT  according  to  (110)  in  order  to  determine 
the  normalized  quality  ratio. 


To  reduce  computation  time,  we  take  advantage  of  various  properties  of  the 
functions  involved.  First,  since  $3(t,  ji)  is  even  in  r  and  (see  (108)),  and  w,  is  even, 
we  express  (F-5)  as 


where 

1/2  for  k  =  0 
1  for  k>l 

Also,  from  (F-l)  and  (F-2),  we  have  normalized  form 


(F-6) 


(F-7) 


'i*3 (Li  a>  Lib)  =  |vo 11 


(1  -  a  -  b)  +  (-l)k_1  V,  sin [kir(a  +  b)] 

k=l  K 


for  0  <_  a,  0  <,  b,  a  +  b  <_  1 


(F-8) 


where  now 


C  =  cos(Tra),  C  =  cos(iTb)  , 
ot  3 

(F-9) 


and  { Vk}$are  still  as  given  in  (F-2). 

Since  L,^3(L,x,  pqL,)  in  (F-6)  is  zero  if  x  +  qk  >  1,  we  can  limit  the  sum  on  k  in 
(F-6)  to  km  =  min  (P  -  1,  1/q),  and  we  need  to  evaluate  the  integral  on  x  in  (F-6) 
only  up  to  Xm  =  min  (L2/L|(  1  -  qk).  The  number  of  x  intervals  needed  in  (F-6)>is 
about  16  with  the  Trapezoidal  rule  for  integration.  These  features  are  incorporated 
in  the  program  listing  below. 
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1  *  NUSC  TECHNICAL  REPORT  6459,  1*9  MAY  1 99 1  ,  H.  H.  NUTTALL 

10  let =10©  !  Be  T 

20  C  1  =  1  ■  6  +  35,  <  48*PI  -  !  HANNING  7Ei<*0FAL-UE  I GHT  I NG  hi  ;  TABLE  I 

30  Q 1 =  3  8  !  62.5*.  OVERLAP  FOR  U1  5  0-01=  =1 

4 0  P fn  =  1  +  v  2 * C  1  * E e t  -  1  1  Q 1  !  P max «  EC'.  1  4  i 

50  PRINT  "Be  T  *,,;Bet,"q  =  q  *  w  1  >  =  "  ;Ql 

60  PRINT 

70  COM  T1,T2,T3,T4 

80  T  2  *  2  *  P I 

90  T3  =  1  • '3 

1O0  T4=l-T2 

110  FOR  P=1  T 0  P m 

120  Be  1  1  =  E e  t  <.  1  +  ( P **  1  ) * Q 1  •*  !  Be  L  1  \  E U •  1 4 U 

130  L21 1=FNL21 1 (Bell , Cl >  !  L2-L1;  SOLUTION  OF  EQ.  52 

140  T 5  =  .  5 * P  I  L 2 1  1 

150  K m  =  M  I N  (  P •  1  ,  101) 

160  S * . 5  *  F  N I n i < O , L  2 1 1 , T5 > 

170  FOR  K=1  TO  Km 

180  S*S  +  <  l-K>P)*FNInt  <Gl*Kf L21  1  ,  T5> 

190  NEXT  K 

200  Q = 4 / P * S 

210  PRINT  P,L2U,Q*Bet 

220  NEXT  P 

230  END 

240  ! 

250  DEF  FNL21 1 (Bell , Cl >  !  SOLVE  EQ.  52  FOR  L2-L1 

260  Eps=E*1 1 - . 5  ' C 1 
270  X 1  =  4  x  C3  *  £  e  1  1  ) 

280  IF  Eel  1  \  1 0  THEN  XI  = .  603245.' SQR <C Eps > -Eps ./ < 6+ 1 2*Eps  > 

290  X2*X 1*1 . 037 
3O0  F 1 =  F  N  F  <  X 1 , Bell  > 

310  F  2 = F  N  F  <  X  2 , E  $  1 1 > 

320  IF  ABS<F2-F1>< IE-6  THEN  330 
330  T=X2 

340  X2*<F2*X1  -F 1  *X2  >  ''<F2-F  1  > 

350  X 1 =T 

360  F 1 =  F  2 
370  GOTO  310 
380  RETURN  X2 
390  FNENZH 
400  ! 
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410  DEF  FNF>: L21  1 , Be  11  >  !  RIGHT  SIDE  -  LEFT  SIDE  OF  EO.  F 

420  COM  T 1 

430  X tit  =  M  I  N  1 , L21  1  > 

440  Tl=. 5*PI ■ L  2 1  1 
4S0  N=16 

4  SO  Del  =  X  fn  N 

470  F=.  5 

480  FOR  K*1  TO  N-l 
490  F  «F  +  FNG ( K *  B  e  1 ) 

500  NEXT  K 

510  RETURN  2*B« 1 *F- 1 • Be  1  1 
520  FNEND 
530  ! 

540  DEF  F N G < X >  !  t  M2  <  L 1  * X > * P h i  1  f  L  1  * X >  1  •'  2 

550  COM  T1,T2,T3,T4 

560  P  *  T  2  *  X 

570  G  a < 1 - X > *  < 2  +  C  0 8 < P  )  )  *  T  3  +  T  4  *  S  I  N  < P  > 

580  W2*C0S< T 1 *X >  "2  !  HANNING  lAG-MEIGHT I NG  M2 

590  RETURN  < M2*G > ' 2 
600  FNEND 
610  ! 

620  DEF  FN I nt  < Qk , L2 1  1 , T5  >  !  INTEGRAL  OF  EO.  F-5 

630  X iri « M I N < L 2 1  1  ,  1-Qk  > 

64  0  Be  I  =  >■; rri -•**  1  6 

650  S«.5*FNPhi3<0, Qk  > 

660  FOR  J=1  TO  15 
678  X*Be 1  *  J 

680  S»S  +  CG$ <  T  5  *  X >  4*FNPh  i  3 < X,  Qk  :• 

690  NEXT  J 
700  RETURN  Bel*S 
710  FNEND 
720  ' 

730  DEF  FMPh i  3<  A, B  >  !  L 1 *PH 1 3 C LI *A , L 1*£ >  for  A>  =  0,  £>«0: 
740  IF  A  +  B  > *  1  THEN  RETURN  0 
750  Ab=P I *< A+B> 

760  C  a=  C  0  S  <  P  I  *  A  ) 

778  C  b  =  C 0 S  <  P I  * B 1 

780  Ca2  =  C  a "'2 

790  Cb2=Cb‘ 2 
800  S*Ca2+Cb2 

3 1 0  V  0  = .  375-8+4*  C  *2  *  C  b  2 + C  a2  -  2  +  C  b  2  2 

820  V 1 *C a*Cb* < 4*8- 1 > 

830  V2=-. 25+ . 5*8+C*2*Cb2 

840  V3«Ca*Cb^3 

850  8  =  VI  *SINt  Ab  >-V2*S I  N».  2*Ab /•  +  V3*S INC  3*Ab  - .  03 1 2 5* 8  In-  4 -Ad 

8  60  L  1  ph  i  3  *  4  <  9  *  P  I  >  *  ^  V  0  *  <  P I  -  A  b  >  +  8  > 

870  RETURN  LI  phi  3 
S80  FNEND 
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